{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Global Forecasting System - Temperature profiles in the atmosphere"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from datetime import datetime\n",
    "import xarray\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from mikeio import Dfs3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's try to download todays forecast from the OpenDAP server."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "now = datetime.now()\n",
    "\n",
    "forecast = datetime(now.year,now.month,now.day)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "dtstr = forecast.strftime(\"%Y%m%d\")\n",
    "hour = \"00\" # valid options are 00,06,12,18\n",
    "url = f\"https://nomads.ncep.noaa.gov:9090/dods/gfs_0p25/gfs{dtstr}/gfs_0p25_{hour}z\"\n",
    "ds = xarray.open_dataset(url)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create a small subset of the data to make it faster to download."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds = ds.sel(lon=slice(10,15), lat=slice(54,58)).isel(time=slice(0,2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Temperature on pressure levels in the atmosphere is named *tmpprs*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre>&lt;xarray.DataArray &#x27;tmpprs&#x27; (time: 2, lev: 34, lat: 17, lon: 21)&gt;\n",
       "array([[[[285.8049 , ..., 284.40488],\n",
       "         ...,\n",
       "         [287.50488, ..., 280.8049 ]],\n",
       "\n",
       "        ...,\n",
       "\n",
       "        [[252.2    , ..., 253.1    ],\n",
       "         ...,\n",
       "         [252.40001, ..., 252.8    ]]],\n",
       "\n",
       "\n",
       "       [[[285.2049 , ..., 282.8049 ],\n",
       "         ...,\n",
       "         [286.8049 , ..., 280.2049 ]],\n",
       "\n",
       "        ...,\n",
       "\n",
       "        [[252.70523, ..., 252.70523],\n",
       "         ...,\n",
       "         [253.20523, ..., 253.00523]]]], dtype=float32)\n",
       "Coordinates:\n",
       "  * time     (time) datetime64[ns] 2020-06-04 2020-06-04T03:00:00\n",
       "  * lev      (lev) float64 1e+03 975.0 950.0 925.0 900.0 ... 5.0 3.0 2.0 1.0 0.4\n",
       "  * lat      (lat) float64 54.0 54.25 54.5 54.75 55.0 ... 57.25 57.5 57.75 58.0\n",
       "  * lon      (lon) float64 10.0 10.25 10.5 10.75 11.0 ... 14.25 14.5 14.75 15.0\n",
       "Attributes:\n",
       "    long_name:  ** (1000 975 950 925 900.. 5 3 2 1 40) temperature [k] </pre>"
      ],
      "text/plain": [
       "<xarray.DataArray 'tmpprs' (time: 2, lev: 34, lat: 17, lon: 21)>\n",
       "array([[[[285.8049 , ..., 284.40488],\n",
       "         ...,\n",
       "         [287.50488, ..., 280.8049 ]],\n",
       "\n",
       "        ...,\n",
       "\n",
       "        [[252.2    , ..., 253.1    ],\n",
       "         ...,\n",
       "         [252.40001, ..., 252.8    ]]],\n",
       "\n",
       "\n",
       "       [[[285.2049 , ..., 282.8049 ],\n",
       "         ...,\n",
       "         [286.8049 , ..., 280.2049 ]],\n",
       "\n",
       "        ...,\n",
       "\n",
       "        [[252.70523, ..., 252.70523],\n",
       "         ...,\n",
       "         [253.20523, ..., 253.00523]]]], dtype=float32)\n",
       "Coordinates:\n",
       "  * time     (time) datetime64[ns] 2020-06-04 2020-06-04T03:00:00\n",
       "  * lev      (lev) float64 1e+03 975.0 950.0 925.0 900.0 ... 5.0 3.0 2.0 1.0 0.4\n",
       "  * lat      (lat) float64 54.0 54.25 54.5 54.75 55.0 ... 57.25 57.5 57.75 58.0\n",
       "  * lon      (lon) float64 10.0 10.25 10.5 10.75 11.0 ... 14.25 14.5 14.75 15.0\n",
       "Attributes:\n",
       "    long_name:  ** (1000 975 950 925 900.. 5 3 2 1 40) temperature [k] "
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ds.tmpprs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre>&lt;xarray.DataArray &#x27;lev&#x27; (lev: 34)&gt;\n",
       "array([1.00e+03, 9.75e+02, 9.50e+02, 9.25e+02, 9.00e+02, 8.50e+02, 8.00e+02,\n",
       "       7.50e+02, 7.00e+02, 6.50e+02, 6.00e+02, 5.50e+02, 5.00e+02, 4.50e+02,\n",
       "       4.00e+02, 3.50e+02, 3.00e+02, 2.50e+02, 2.00e+02, 1.50e+02, 1.00e+02,\n",
       "       7.00e+01, 5.00e+01, 4.00e+01, 3.00e+01, 2.00e+01, 1.50e+01, 1.00e+01,\n",
       "       7.00e+00, 5.00e+00, 3.00e+00, 2.00e+00, 1.00e+00, 4.00e-01])\n",
       "Coordinates:\n",
       "  * lev      (lev) float64 1e+03 975.0 950.0 925.0 900.0 ... 5.0 3.0 2.0 1.0 0.4\n",
       "Attributes:\n",
       "    grads_dim:      z\n",
       "    grads_mapping:  levels\n",
       "    units:          millibar\n",
       "    long_name:      altitude\n",
       "    minimum:        1000.0\n",
       "    maximum:        0.4\n",
       "    resolution:     30.290909</pre>"
      ],
      "text/plain": [
       "<xarray.DataArray 'lev' (lev: 34)>\n",
       "array([1.00e+03, 9.75e+02, 9.50e+02, 9.25e+02, 9.00e+02, 8.50e+02, 8.00e+02,\n",
       "       7.50e+02, 7.00e+02, 6.50e+02, 6.00e+02, 5.50e+02, 5.00e+02, 4.50e+02,\n",
       "       4.00e+02, 3.50e+02, 3.00e+02, 2.50e+02, 2.00e+02, 1.50e+02, 1.00e+02,\n",
       "       7.00e+01, 5.00e+01, 4.00e+01, 3.00e+01, 2.00e+01, 1.50e+01, 1.00e+01,\n",
       "       7.00e+00, 5.00e+00, 3.00e+00, 2.00e+00, 1.00e+00, 4.00e-01])\n",
       "Coordinates:\n",
       "  * lev      (lev) float64 1e+03 975.0 950.0 925.0 900.0 ... 5.0 3.0 2.0 1.0 0.4\n",
       "Attributes:\n",
       "    grads_dim:      z\n",
       "    grads_mapping:  levels\n",
       "    units:          millibar\n",
       "    long_name:      altitude\n",
       "    minimum:        1000.0\n",
       "    maximum:        0.4\n",
       "    resolution:     30.290909"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ds.lev"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Level 0 = 1000 mbar, i.e. close to the ground."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.QuadMesh at 0x20cc61dfac8>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZUAAAEXCAYAAABlI9noAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOydebwdVZW2nzdzSEIghEHGAI3SoBAw0mocQBGRFgdURFFwxAG6AbUdsFsFh1ZAbdp2igZxAAcMET6USWUQFSGBQMAgqARMwhSSkISM9973+2PvS05OzlB1cu69Jzfrya9+t2rXXrV3nXNSq/Zea68l2wRBEARBOxgy0B0IgiAIBg+hVIIgCIK2EUolCIIgaBuhVIIgCIK2EUolCIIgaBuhVIIgCIK2EUqlH5G0p6SVkoYOdF+C9iBpkiRLGjbQfQmCTiCUSh8iab6kI3uPbT9ke6zt7oHsVz0kjZQ0XdKDklZIukPSq6rqvFzSvZJWSbpe0l4V586XdH+WvVfSSVWykyXNzrKzJU1u0p+6beXzR0q6XdJTkv4h6fgm93ahpOWSHpH0oTr1Pp2VxJG1znc6kp4t6RpJiyVtsghN0gRJM/Nn9qCkt1adf2suf0rSLyRNKCpbdR1J+pKkJ/J2riS1926DTiSUSlDJMOAfwEuB8cB/AT+TNAlA0kTgslw+AZgF/LRC/ing2Cx7MnCBpBdm2RHA5cCPgO2B7wOX5/JNaNaWpAOAS4BP5vYmA7Mb3NtngP2AvYAjgI9KOrqqzX2BNwIPN7hOp7Me+Bnw7jrnvw6sA3YGTgS+KelAgPz328Db8/lVwDeKyNbgFOB1wMHAQcCrgfe1fFfBloPt2PpgA34I9ACrgZXAR4FJgIFhuc4NwOeAP+Q6/w/YAbgYWA7cBkyquOb+wHXAEuAvwPH9cB93AW/I+6cAf6g4Nybf3/51ZK8APpz3jwIWAqo4/xBwdB3Zhm2RFMpnS9zHQuCoiuPPAj+pqnMVcAwwHziy4HWrv9PxwHSSYlqYv9+hwEhgGfDsCtkd8z3t1Aff2z+l/94blY0hKYVnVv1Ov5j3vwBcUnFu31x/XDPZGu3/ATil4vjdwC19/XuNbeC3GKn0EbbfTnpoHus05XVunaonkN4MdyP9J/4j8D3S2/k84NMAksaQFMolwE7AW4Bv1HtTlPQNScvqbHcVuQdJOwPPBO7JRQcCd1bc41PA33J5texo4HlVsnc5P2Eyd9WSLdjW83M7cyU9LOlHlVM1VX3ZHti18np5/8CKOm8C1tn+VZ3+FOX7QBfpoX4ISZm+x/Za0sjrLRV1jwdutP1YjT6/qMH3t0zSi1ro2zOBbtv3VZRVfg7Vn/nfyIqkgGw1G12rSd1gEBFKZeD5nu2/2X6S9Kb8N9u/tt0FXEp6MEGaPphv+3u2u2zfDswgTddsgu0P2t6uznZQs05JGk4aMX3f9r25eCzwZFXVJ0lvstV8i/QguaYF2SL1dycp4zeQprVGA19rcK1e+U2uJWks6S39jDryhchK+FXAGbafysriq6QXB0gvBJVK5a25bBNs39zg+9vO9s0tdLHZZ9ro/OZ+f08CY8OuMvgJj5WB59GK/dU1jnsfiHsB/yJpWcX5YaQpiLYiaUi+7jrgtIpTK4Ftq6pvC6yokj8PeDZwRMXIpK6spD2BP/cW2h5boK3VJIV8X27zC8Cv8/63gLflel8Avlkhv6bGtc4Gfmj7ATaPvYDhwMMVz84hJDsVwG+B0ZL+BXiEZAeauZltlqHZZ9rofE8T2WZtbQusrBqpBoOQGKn0Le38D/QP0lRJ5dvqWNsfqFVZ0reU3JdrbffUkslyItkEdibZUtZXnL6HZHjtrTuGNGV3T0XZ2aS39aNsL6+SPajqTfUg4B5v8IobmxVKkbbuos7na/v9Fdf7gu2lJBvHwRXVDq641suBf89eYY8Ae5AcFD5W52Oqxz+AtcDEiu9oW9sH5n71kIzobyGNUq60XfOhLOnFDb6/lZJeXLJvAPcBwyTtV1FW+TlUf+b7kGxB9xWQrWajazWpGwwmBtqoM5g34BY2NlZOYlND/Xsqzn8OuKji+Ejgr3l/HPAgacpneN6eB/xzm/v8rdzvsTXO7UiaxngDMAr4EhXGV+ATwP3AM2rIjsj9P530oDotH4+o049mbb0LeADYB9iG9LD+YYP7+iJwI8nzbH+Skjk6n9sB2KVi+wfwpt7PgOQ5dkOd61Z/p5cDF5DezIeQFOFLK+r/S277buC1ffCbU/68Dsj9GgWMrDj/E+DHJMP71PwZH5jPHUhyEHlxPv8jKpwZGsnW6Mf7STbB3Uj2rHuA9w/0/8nY+n4b8A4M5g14LclYvwz4SI0H0A0UVCr5+FnAL4HHgSdI0ymT29jfvXL/1pCmL3q3E6v6dC9p+ukGNvZOM+lNvVL2rIrzh5DcflcDtwOHNOlP3bby+bPzZ/E4abpu+wbXGglcmB+ajwIfalB3PhXeX6SR2+fr1K3+TseTptsW5IfuHcAJVTJ/JXnw1VSom/kd9vancptfcX4C8AuS+/dDwFur5N+ay58iKcgJRWRJimhlxbGAc/N9Lsn7avf9xtZ5m/IPIAiCOkiaA7zc9hMD3Zcg6HRCqQRBEARtIwz1QRAEQdsIpRIEQRC0jVAqQRAEQdsYVIsft5sw1Lvu3rdR5Ye1sPRkaAuLiEV5mSEtyHTTU6r+492jS7exePXY5pWqsFv4zIaU+27UwnfZynrw7u4W3t3KfS2t0d3CzQxt4fe/ulw7w1d0lW6DrvIf2PL1jy22vWP5xjbwyiPG+IklxYKOz75r7TW2j653XtIewA9Ibu09wDTbFyhF8/4WyT28C/ig7VsljSe5fe9Jepafb/t7m3M/7WBQKZVddx/Kj67cpU/bmDBkffNKVWw3pPzHPLKF9BwjNby0zPKe1aXqT1v27NJtfOeeqaVlutaVfzkYMbrcdzNsaPkH0bCh5bMWLH9ym9IyLnn/bkFBDF1Z/jPuGV/+gb/tXeV+l7v+pryTnZbWW9hfn6sX/O+DpYWqeGJJN7des2ehukOfcf/EJlW6SAFYb5c0Dpgt6TqSO/bZtq+SdEw+Phw4Ffiz7WMl7Qj8RdLFtte1ej/tYFAplSAIgv7EQE+bhpW2HyanXbC9QlLv4lGzIeTNeGBRRfPjcpSKsaT1QC0M89pLKJUgCIIWMWZ98Zx7EyXNqjieZntarYo5h9EhwJ9IgU6vkXQ+yQ7+wlzt/0jpJRaRIm682SkU0IDS54Z6peyHcyXN6f1AlTIA3tJbJumwOrInK2USvF/SyX3d1yAIgrL0FPwHLLY9pWKrp1DGkiKQn+EUP+8DwJm29wDOJEV4AHglMIcUBmcy8H+SqoN+9jv95f11hO3Jtqfk4945wsnAp/LxRuTcGJ8mxUo6DPh0zosRBEHQERjT7WJbEXLKiRnAxbYvy8Unk3LxQEqH0fsS/k7gMif+SoqFt3/bbq5FBsqluN4cYSWvBK6zvcQpyux1QF3PiSAIgoGgBxfamlERIXye7a9UnFpESvEN8DJS0FZI8ddenmV3JsUG/Hubbqtl+sOmYuBaSQa+nYd89eYIK9mNDXkoIAXo2626kqRTSKln2WW3vnUnDoIgqMRAd/syXEwlRSGfm+PNAZwFvBe4QNIwUrDXU/K5zwIXSZpLCuD5MduL29WZVukPpTLV9iJJOwHXSbqXlK3wTNszJB1P0s5HVsnV8pHc5NvLSmoawAEHjYhAZkEQ9CtFRiFFcMrmWc83/Lk16i8ipavuKPp8+ivfOE6pVWeS5gPrzRFWsoCULKmX3ak9TRYEQTAgGFhvF9q2FvpUqUgakxfx9GbuO4qUnKjeHGEl1wBHSdo+G+iPYkO+8yAIggHHmO6C29ZCX09/7QzMzBlkhwGX2L5a0kpqzBFKmkLKDvce20skfRa4LV/rHNtL+ri/QRAExTF0bz36ohB9qlRs/52N81T3lt9M7TnCWcB7Ko4vJGXrC4Ig6DjSivqgklhRHwRB0DKiu4VAroOZQaVUVvWM5PY1e5WSGalyQQgPHFneV2C7ktFzAda6fAifVmTKBq6csk15N/grtn9OaZlWeHx5uWjITz1ZPuIyq8u7rWts+SCkI8evKVV//UPlI0EPX1H+YdjVXf6Rsc1j5d7lu8eOKt3GsBaiFLcDAz0x/bURg0qpBEEQ9CcG1kVaqo0IpRIEQbAZ9LSQ+2cwE0olCIKgRdKK+lAqlYRSCYIgaBEjumP6ayNCqQRBEGwGMf21MaFUgiAIWsSIdY5AtpWEUgmCIGiRtPgxpr8qCaUSBEGwGYShfmNCqQRBELSILbodI5VKQqkEQRBsBj0xUtmIUCpBEAQtktapxEilklAqQRAELWLEesdjtJJB9Wks7dqGGY8cWkpm7LB15RrZsVx1gF2HPlReqJ9Y1lMuCOU+w8oFOgQ4aEL5IJxXzTuwtIxXDi9Vf8jq8m+YPcPLRw/sWd3Cf7O/lwuqOO7x8k0MLfnTB+gaXX6qZ/025WR6RpZ30fWwgRstdMc6lY2IcVsQBEGL9K6oL7I1Q9Iekq6XNE/SPZJOz+WTJd0iaY6kWZIOy+X/kcvmSLpbUrekCX18y03p85GKpPnACqAb6LI9RdJPgWflKtsBy2xPLiLb1/0NgiAoQ0/7vL+6gA/bvj2nYZ8t6TrgXOBs21dJOiYfH277POA8AEnHAmd2Qnbc/pr+OsL24t4D22/u3Zf0ZeDJorJBEASdQjsN9bYfBh7O+yskzQN2y81sm6uNB2rNJ78F+HFbOrKZDKhNRSl5/fHAywayH0EQBK1gVMamMlHSrIrjaban1aooaRJwCPAn4AzgGknnk0wWL6yquw1wNHBaqc73Ef2hVAxcK8nAt6s+xBcDj9q+vwVZACSdApwCMGrnce3teRAEQQNsynh/LS4yhS9pLDADOMP2ckmfI01tzZB0PDAdOLJC5Fjg950w9QX9o1Sm2l4kaSfgOkn32r4pn2s2ZGskC0BWNNMAxu+/cyT2DIKgH1FbFz9KGk5SKBfbviwXnwycnvcvBb5bJXYCDZ6jknYBsP2IpB1JL/N/sX1P2zpeQZ97f9lelP8+BswEej0XhgHHAT8tKxsEQdAJGOj2kEJbM7I5YDowz/ZXKk4tAl6a918G3F8hMz6fu7zONd8H/BG4RdIHgCuBVwOXSXp36RsuQJ+OVCSNAYZko9MY4CjgnHz6SOBe2wtakA2CIOgI2riifirwdmCupDm57CzgvcAF+UV8DXm6P/N64FrbT9W55mnAgcBo4EHgn/KIZXvgepISayt9Pf21MzAzKWCGAZfYvjqf22TIJmlX4Lu2j2kiGwRBMOAYtS1Jl+2boe5c2nPryFwEXNTgsuttrwJWSfqb7Uey3NJsq247fapUbP8dOLjOuXfUKFsEHNNMNgiCoBMwpQz1A0GPpOG21wP/2lsoaRR9ZP7o6E8jCIKgs1Gn51M5jqT7qDI17AB8uC8aHFRKpbtnCEtXb1NKZuy4cgGQHli7U6n6AEtG/aO0TH+xpuRb1sKu7Uq3scuI5aVl3F3+JWrkY+ViRnWPKN0Ew1aVf4AMWV/+XoauKlm/hTheQ9eWn/0YurZ8O0O6y9X3sPKfsUeWi/vWLkxbV9S3Hds1Aw/aXggs7Is2B5VSCYIg6G86fKTS74RSCYIgaBFbHT1SGQhCqQRBEGwGkU54Y+LTCIIgaJGUpGtooW0gkHR0xf54SdMl3SXpEkk790WboVSCIAhaJBnqVWgbIL5Qsf9lUhTkY4HbgG/3RYMx/RUEQbAZbEE56qdU5K36qqST+6KRUCpBEAQt0s4V9X3ETpI+RFqpv60k2e71JY/Fj0EQBJ1GT2ePVL4D9OYE+T4wEXg8Ry6eU1dqMwilEgRB0CI2ZZJ09Tu2z65T/ghwUl+0GUolCIKgRYzo6hkYz65OJZRKEATBZhAr6jemoVKR1Cxok4CHbT+zfV0KgiDYMuh1KQ420Gyk8jfbhzSqIOmONvZnsxiiHsaMKBfxbsdRK0rVP2Sb+aXqAyzrGVlaphXWuHxQvSe6x5aqv2j99qXbGDVkfWmZg/cpH4Tzr/fsW6r+mCdLN8HwVeWDMHaNLP/QWTe+XP2lB/aUbmPoqvIG5ol3lW9nzKJy0S6HLVtTug11l4xa2Ta23DAtknbpza9SUXZoAdH1tufWO9lMqbyhQANF6gRBEAxK2pmjvp+ZTkWOlcyNpIWRjW5qb2BSvZMNlUpOlNWQInWCIAgGIzas30IN9barFQrAbbZf1khO0m8bnS80bpN0nKT7JT0pabmkFQXsLb2y8yXNlTRH0qxc9tN8PCefr+kvLeloSX+R9FdJHy/SXhAEQX/Ru/ixHWFaJO0h6XpJ8yTdI+n0XD5Z0i29z1BJh1XIHJ7L75F042bfTxOFUqROUe+vc4Fjbc8rWL+aI2wvrujUm3v3JX0Z2GR2W9JQ4OvAK4AFwG2SrrD95xb7EARB0HbaOP3VBXzY9u2SxgGzJV1Hev6ebfsqScfk48MlbQd8Azja9kOSSmUQlHSl7VfXOXeO7U9VHA8FfmD7xGbXLWphenQzFEpdJAk4HvhxjdOHAX+1/Xfb64CfAK9tdx+CIAhapZ0BJW0/bPv2vL8CmAfslpvZNlcbDyzK+28FLuvN7mj7sZLdf2+Dc3tK+gSApJHATOD+Ihdt5lJ8XN6dJemnwC+Ap92rbF9WoA0D10oy8G3b0yrOvZiksGp1djeg0gVoAfAvNfp4CnAKwMidx1WfDoIg6FNKeH9N7DUBZKZVPQ+fRtIk4BDgT8AZwDWSzicNBF6Yqz0TGC7pBlIolgts/6Be45ImALa9lLTzcIO+vhO4OCuWI4CrbH+16R3SfPrr2Ir9VcBRFccGiiiVqbYX5aHZdZLutX1TPvcWao9SoLb3wSb+nPlLmQaw7bN2Lu/vGQRB0Crlwtovtj2lWSVJY4EZwBm2l0v6HHCm7RmSjid5bR1Jen4/F3g5MBr4o6RbbN9Xca09SdNlLweWpSJtC/wW+Ljt+VVtV7oUX0AKj/974EZJh/aOpBrRzPvrnbmhqbZ/X9X41GYXz9dYlP8+JmkmaVrrJknDgONIH0otFgB7VBzvzoZhXxAEwYBjoKuN61QkDScplIsrZoJOBk7P+5cC3837C0iK6ingKUk3AQcD91Vc8qfA/wAn2u7ObQwF3kQyKTy/qgtfrjpeChyQyw00NeQXNdR/DaheFFOrbCMkjQGG2F6R948CzsmnjwTutb2gjvhtwH6S9gYWAieQ5hCDIAg6gnauqM825unAPNtfqTi1CHgpcAPpod5rLrgc+L/8gj6CZB6onqKaaPunG/U5KZefSPpsdR9sH7G599HMpvIC0vzdjjkmfy/bAkWcs3cGZqbPimHAJbavzudOoGrqS9KuwHdtH2O7S9JpwDW5rQtt31OgzSAIgn6jjWFapgJvB+ZWLLM4i2RQvyArjzVkG7LteZKuBu4CekjPzrurrjlb0jdIYe97bdR7kEY/m0RDkfRq21c26mSzOs1GKiOAsblepRV8OfDGJrK9CyMPrnPuHTXKFgHHVBz/CvhVs3aCIAgGgnYm6bJ9M/VXstc0E9g+DzivwWVPAt4NnE1yfhJp2uwK0qiomvMkLWzQD0gpiltTKrZvlHQz8Jx6cfmDIAi2Zjo5TEtejvHNvBXhUeArTeo0dC1ualOx3Z1d0TqetWuH87f5u5QTmlSu+gvGbVNOgNYCPe4wdGVpmYUtBHtcUjKg5Nqe8veypgWZVhj24iWl6q+8q/zPesI95R0Mt3m8fLDDcQvKtbNDCxPDw54qH7hx2LLVpWW0ulxASdaUCwo7oLizoxTnKbN3A69jw5qXRSR7zHTbG0V7tX345rZZ1FB/h6QrSJ4HT1V0oIhLcRAEwaDEQFdPR0cp/iHJlfhs0rQXJE/ak4EfAW+uI9cyRZXKBOAJNnYnK7pOJQiCYFDSTptKH3Go7WdVlS0AbpF0Xy2BzaWQUuldrxIEQRBsjDtbqSyV9CZghu0eAElDSOtUlvZFg0WjFO8uaaakxyQ9KmmGpN37okNBEARbEj2o0DZAnEDy1H1U0n2S7gceIS08P6GekKRtJP2XpO/k4/0k1Qw+WU3R6a/vAZeQtBvA23LZKwrKB0EQDDrc4Yb6HIblzQCSdgBUGTG+Ad8DZgMvyMcLSDb1hmtYoLhS2dH29yqOL5J0RkHZIAiCQYro7mxDfW8ssaNJix678mjl2t7psDrsa/vNkt4CYHt1XvHflKKfxmJJb5M0NG9vIxnugyAItmpsFdoGghyA8nqSUjmNFHvx7cAcSc9pILpO0mhyEF9J+1IRob4RRUcq7wL+jxRXxsAfclkQBMFWSztjf/UR/wk83/YqSRNJgSpfKekgUgTiF9aR+zRwNbCHpItJIWTeUaTBot5fDwGvKVI3CIJgq8HJrtLBCOhdsfoUsBOA7btyCPxNBdI0170kY/7z8zVOL2iLKaZUJO1ICmo2qVLGdoxWgiDYqunkMC2k2IlX5/z1ryIZ23sTdtXsuG1L+oXt5wK/LNtg0emvy4HfAb8GysecCIIgGISYzl6nYvtjOa/9AcA5tq/Lp5bROHXJLZKeZ/u2sm0WVSrb2P5Y2YsHQRAMbkR3T+cqFagd7T17fjUyvB8BvE/Sg6RpMyUxH9SsvaJK5UpJx+TOdS49glVF0rxsYMHS7UrV//nwptlAN+GJNeWDUL5vr5uaV6pibQuBK5/sKte3pSXrAyxeWy5oJcCjT41rXqmKlStHl6o/snzMTrabW34RsoeW+00CDFm7vnmlzWV9V3mZFoI9evmKcvW7WujXANLJI5XN4FWtChZVKqcDZ0laC6xng9aqaegJgiDYGrAHrVJp2f2gqPdXw9dGSQfWy8ooaT6wgmSL6bI9JZf/G8lvugv4pe2PFpUNgiDoFDrcpbhVfklSLAJGAXsDfwEObCZYdKTSjB/S2OhzRKU7mqQjgNcCB9leK2mnorJBEASdRCe7FEs6ujeFu6TxpARczwPuBs60/WgtOdsbLYyUdCjwviJttiu+QFlV/QHgi7bXAth+rE39CIIg6DeM6OkZUmgbIL5Qsf9l4GHgWOA20uLHQti+naSMmtKukUojXW3gWkkGvm17GvBM4MWSPg+sAT5Sx3WtlmwQBEHH0K6BiqQ9gB8AuwA9wDTbF0iaDHyLNA3VBXzQ9q2SDict93ggX+Iy2+c0aGKK7cl5/6uSTm7Qlw9VHA4hzUQ9XuQ+2qVUGjHV9qI8xXWdpHtzu9uTVms+D/iZpH3sTQaSm8ja3sgtStIpwCkAQ7cvn043CIKgZdprqO8CPmz7dknjgNmSrgPOBc62fVVec3IucHiW+Z3tRiHpd8oKQsC2klTxnG00fKq0o3eRbCwzitxEu5RK3STUthflv49JmkkKaLaApFUN3CqpB5hIlSasI3tTVZ1pwDSAkXvu0cGzm0EQDEra9NSx/TBpegrbKyTNY0Ne+V5P2/GkHPNF+Q4bFMT3yc9ZSbsAcxrI/dn2pZUFOdnXpXXqP03RMC1TgTm2n8oRig8FLrD9IIDt59eRGwMMyR/QGOAo4BxgJSk18Q2SngmMABYXlA2CIOgYSoxUJkqaVXE8rd6UvqRJwCHAn4AzgGsknU8aXVQGgXyBpDtJiuYj1V64ts+u3Wc/ApzUoK+fYFMFUqtsE4qOVL4JHCzpYOCjwHTS3N9Lm8jtDMzMYfiHAZfYvlrSCOBCSXeTRjkn53gzuwLftX1MPdmC/Q2CIOgXSnh/LS6yLCLnP5kBnGF7uaTPkTy1ZuRQ9tOBI4Hbgb1sr8zTYr8A9qtxvf1JI54/2V5ZUX509TNV0quAY4DdJP1vxaltSdNgTSmqVLryQ/+1pBHK9EZGnl5s/x04uEb5OlL2yOryRaQbqisbBEHQKdjgNnp2SRpOUigX274sF59MWoAOaaTw3dS2l2/oh38l6RuSJlYt3/h34FRgHjBd0um2L8+nv0AKb1/JImAWKSr97IryFcCZRe6hqFJZIekTpOQuL5Y0FCgfEyQIgmCQ0a51Kjnk/HRgnu2vVJxaRJoVuoFkNrg/198FeDS/8B9GmhqrTp74XuC5eTQzCfi5pEm2L6DGUhDbdwJ3SrrEdkuxgooqlTcDbwXeZfsRSXsC57XSYF8yZEQ3Y3YvF2fogB1rrv2py7qe8nGcdhi1qrTMku7y8bJaYW1POV+NBavKe9j9Y0W5+GoAjy8tH/uLhaNKVd/pjn6IrwUMWbm6eaVqVj5VqrpXtdBGC7QSl8slY4y5q/z3MnRcC7+XdtE+96CppBf3uZJ6jehnkRTDBZKGkZZgnJLPvRH4gKQuUs6UE2p40A7tnfKyPT+7If9c0l40Xl84SdJ/k6IbP/0fy/Y+zW6iaJiWRyTNYMN83WJgZhHZIAiCwUv7UgXbvpn6D/rn1qj/f6SMvI14RNJk23OyzEpJrwYuBBqlE/4eKfvjV0kRi9/ZoG8bUWgyUNJ7gZ+zYQXmbiSjUBAEwdaNC24Dw0nAI5UFtrtsnwS8pIHcaNu/AWT7QdufIU29NaXo3MeppDUif8qdur9JvK4gCILBT4dHKba9AJ7O3rs7yYPrAdsrbf++gegaSUOA+yWdBiwkpyJuRlG3hbXZY4vcwWEMpO4NgiDoFKxi2wAg6QBJvwb+SBoUfJdks7koB5isxxnANsC/k6be3kbyQmtKUaVyo6SzgNGSXkFya/t/BWWDIAgGL509/XUhcKrtfwJeBNxre2/g9yRPs03I3r3H59HMAtvvtP0G27cUabCoUvk4KYTKXFL4418B/1lQNgiCYPDS2UpltO2/ANi+lWyct/0dkmfXJtjuBp6bXZxLU9T7q0fSj4CbejsYBEGw1WMGbGqrIH+T9F/Ab4DjyPG+8iLLRs//O4DLJV1KylEPQMWCzLoU9f56Te5Mb7KXyZKuKCIbBEEwmEkphZtvA8S7SAElzwLWsmFl/jY0jv01gbSQ8mWk/CvHAo2iIT9NUe+vT5O8v24AsD0nr84MgiDYuunp3JGK7WWkeI3V5U8CdW0ktt/ZaptFbSpduRNBEARBBXKxbUtC0jMl/SYH/UXSQZIK2dGLKpW7JbNCnlgAACAASURBVL0VGCppP0lfA/7QYn+DIAgGB0WN9FuYUiHlYfkEsB7A9l3ACUUEiyqVfwMOJM3JXQI8SfJjDoIg2IopuEals435tdgme4tV0p7Q99ln+Wzb/wF8soXO9Ru2WLeuXIDE8SPKBeIbN2xNqfpQPmgjwNyVu5eW6VSeWD6mtMx2vxldWmbVLuXqP7lP+UDboxaWFoHhLSRY3b7RurQa7Fg+0GfXduU/48cnb1NaZs2O5eq7fMxW1u5UPtAl7y0vUpMOHoXkhervBl4P7Erq7SJSbvvpDSIRL5a0b66PpDeSs1I2o+mv3Xa3pE2CmQVBEAR0tFIBfggsAz5DSuMOKVzLycCPSBHoa3EqKU37/pIWAg8AJxZpsOgr1B3Zhbi0z3IQBMGgxXS09xdwqO1nVZUtAG6RdF89oZwk8cjKtO5FGyxqU2nZZ1nSfElzJc2pzM8s6d8k/UXSPZLOrSN7dK7zV0kfL9jXIAiCfqPDvb+WSnpTDg6Z+isNkfRmYGk9IUk75HTCvwNukHSBpB2KNFh0RX3LPsuZI6pSXB4BvBY4yPbaWhGPsy3n68ArSJr1NklX2P7zZvYlCIKgfXT29NcJwJeAb0haSsqJMh64nsbeXD8BbgLekI9PBH4KHNmswUJKJWusap4EZlXkOy7DB4Av2l4LYPuxGnUOA/6ah2FI+glJEYVSCYIgKIDt+WS7SR5pqPIFvwETbH+24vhzkl5XpM2i01+jgMmk3Mj3AweRpsTeLel/msgauFbSbEm9aTCfScp1/ydJN0p6Xg253YB/VBwvyGUbIekUSbMkzepeXi4FaxAEwebS4dNfT2P7CWCcpOMk7d+k+vWSTshTZUMkHQ/8skg7RZXKPwEvs/01218jDYH+meSmdlQT2am2DwVeBZwq6SWkEdL2wPOB/wB+ViMiZi3r1yZfje1ptqfYnjJ02/Kuq0EQBJtFm9apSNpD0vWS5mVb8+m5fLKkW3rt0pIOq5J7nqTu7PZbfc1fVOy/FvgtySZ+haR3NOjO+0hrEtfl7SfAhyStkLS80X0U9f7aDRhDmvIi7++a3Y3XNhK0vSj/fUzSTNK01gLgMtsGbpXUA0wkhdfvZQGwR8Xx7iT/6iAIgs7AQE/brtYFfNj27ZLGAbMlXQecS1oreJWkY/Lx4fC07flLwDV1rrlXxf7HSIODByRNJEUuvqiWkO1xrd5EUaVyLjBH0g2kEcRLgC9kd7Nf1xOqdEfL+0cB5wArSZ5kN0h6JjACqJ7nuw3YT9LepFSWJwBvLXpjQRAE/UG7prZsP0xeYJifmfNIL/QGts3VxrPxy/W/ATOAWiYE2Hh2Z5jtB/L1F+eX+bpIOgiYRIWeKLKMpKj313RJvyKNMgSc1TsCIU1f1WNnYGae2RoGXGL7akkjgAtzsLJ1wMm2LWlX4Lu2j7HdlXMjXwMMBS60fU+R/gZBEPQbxZXKxMplFcA029NqVcxR4A8hpQA+A7hG0vkkk8ULc53dSCaIl1FfqRycp6sEjJS0i+1H8jO4buwCSReSbOf3sGEsZqA9SiXbO14O7GP7HEl7SjqsRmyYjcieWwfXKF9HynlcXb4IOKbi+FekLJNBEASdSXGlstj2lGaVJI0ljT7OsL1c0ueAM23PyAbz6SS79v8AH8tmiNpdc92gN9uQ7Cb1eL7tmpkhm1HUUP8N4AXAW/LxCtIakiAIgq2Wop5fRafIckbGGcDFFVNNJ7NhhHApacYIYArwE0nzgTeS1qLUdfuVtH221WB7me0/NujKHyW1pFSK2lT+xfahku7IHVqah08dxajh6zngGY+Uktlz1JJS9Zd2lQ+ot2L9qNIyi1Zt27xSFX/78yYe103Z/s9F3ysSqzdZptqcIeWaAGDlnuVlevYv51K+dEn572XiHeVlhj1c7jcG0PWMCaXqd48pH7TyyUkjy8tMWVdaZsz4ckFbu7rL/2BGtCDTNtoUpiXPCE0H5tn+SsWpRcBLSUkSX0Za1oHtvStkLwKutP2LCjmySeGLpDV+Y4GFeVRzIfD5BgElv09SLI+QotMrNemDmt1H0V/i+uxl0Buxckfa6fMQBEGwhdLGNShTgbcDcyXNyWVnkeIpX6AUcXgNcEod+Vr8CDjH9kmSjgNeDPwnKVfK1xtc68LevlDyWV9UqfwvMBPYSdLnSUOtQlnAgiAIBjXt8/66mdrr8wAaRoq3/Y46p3awfUOuc5mkT9p+CvhPSfc2uORDtq9o0uWaFPX+uljSbJKxXsDrbM9rpcEgCIJBQ4eslm/A45LeRlr0+AZgPjw91dZozvBeSZcA/480/QW0waVYUuXE7mPAjyvP2S4/WRwEQTCY6Gyl8i7gfODjwBzgtFw+gTQFVo/RJGVSGTGlLS7Fs/OFBOxJCpUsYDvgIWDv+qJBEARbAR2sVGw/BBxfo/wJkpdZPbmWI9M3dJmwvbftfUgLEI+1PdH2DqRcKpGgKwiCrZ4tJaBkGSQ9U9Jv8gJ1JB0kqZAdvagf3vPyQkQAbF9FcnELgiDYunHBbcviO6TpsfUAtu+icf6Vpynq/bU4a6kfkT6et5EyQQZBEGy9bIGjkIJsY/vWqpX6XUUEi45U3gLsSHIrnpn339JQIgiCYGtgCxqpSHqRpA9JapayZLGkfdmwNvGN5GCXzSjqUrwEOL1I3SAIgq2KDlEYtZB0q+3D8v57gVNJA4NPSzrU9hfriJ4KTAP2l7QQeICUUrgpDUcqkj5ToNNN6wRBEAxGRMcb6odX7J8CvML22SRX4UZKwraPJM1K7W/7RRSc2Wo2UnlPkyxfIhlvPlOksSAIgkGFoXFWkgFniKTtSQpBth8HsP2UpEY2khnAoXn1fS8/p8nKfmiuVL4DNMsA9p1mjfQXq9eMYM69ezWvWMGfx+9Sqv5L9/5rqfoAE0euLC2zeG351Mgjl5QPqqfukm204J7RNbq8zLop5YJDAkzebWGp+v/YdrvSbTz46vIRNccsLP9d9tTNdFGbugHOG7B+bHkZVpdvaNiEkj+yFlg3kAElO3j6i5TUazY5IGRFPpWx1AgJk3PXHwiMz7HCetkWKBRNtaFSycOkIAiCoB4drFRsT6pzqoeU4KuaZ5HWIW5HymXfywpSYMumlI+XHQRBEDxNJ7sUS9rO9rLqcturSMb36vLLgcslvaBJvpW69PmYUdJ8SXMlzelNpSnpM5IW5rI5ko4pKhsEQdBRdLZL8WJJv5b0bkmF53tbVSjQfyOVI2wvrir7qu3zW5QNgiAYeDrfUD+PlHb4LcC5km4mBQa+3Ha57GkFKTRS2Zw4MEEQBIOazh6prLd9pe0Tgd2Bi0kBJhfk0PZtp+j0V8txYEgf57WSZkuqzDJ2mqS7JF2YXd7KyD6NpFMkzZI0q3tFeY+hIAiCzaHD16k87eFle7Xtn9k+DugNFFxbSJoi6UxJ50k6R9LxValQ6lJUqWxj+9aqskJxYICptg8FXgWcKuklwDeBfYHJpKX/Xy4huxG2p9meYnvK0HHlXTeDIAg2i84eqVxcq9D2k7a/X10u6R2SbicNIkYDfyHl0noRcJ2k70vas1GDZQJKthQHxvai/PcxSTOBw2zfVHET3wGuLCoL3FSrbhAEQb/TQXG9alHQbl3JGNLLfE17i6TJwH6kfFo1KTpSORX4NhviwJwBfKCZkKQxksb17pNCA9wt6RkV1V4P3F1UtmB/gyAI+hyV2JpeS9pD0vWS5km6R9LpuXyypFt6vWAl9cbyem02IfSWv6jGNYdIepekX0q6M5sSfiLp8Fp9sP31RgZ823Ns/6bRfRQNKPl34Mj8cB9ie0UROWBnYGYOnzwMuMT21ZJ+mDWeSTmT3wcgaVfgu7aPqSdbsN0gCIJ+oY3eX13Ah23fnl+oZ0u6DjgXONv2VXn5xbnA4cBvgCtsW9JBwM+A/auuOR14EPhv4I3AcuB3wH9Keo7tr210L9Iw4N2kl/1dSc/oRcDlwHTb65vdRLMc9R+qUw6A7a80ks/K6OAa5W+vU38RcEwj2SAIgo6iTdNfth8mmxVsr5A0D9gtt7Btrjae9JDHdmX8pzF1evLcitTAN0u6xfanJN1Eyln/tar6PwSWkeI5LshluwMnk/JpvbnZfTQbqfTG/XoW8Dzginx8LGHbCIIgKKNUJlYt4p5me1qtipImAYcAfyKZG66RdD7JZPHCinqvJ41CdgL+tcal1kva1/bfJB0KrAOwvVaq6ZN2qO1nVZUtAG6RdF+BeywW+0vStbmxFfn4M8ClRRroTzTUjBi/tpRMT8lAdHc8vlup+gCH7Fgu0CHACyZsEkGhKY9NKR8h8Imd6nlzt4/hE9aUlvmnnR8vLbPjqKKzsolRQ5uO5DdhzeGbRLxoysquEaVl1nW3ECGyJE+tG1laZsT68uul16wb3rxSBWtXlu+XnhygiFPl3IUX257SrFIO9jgDOMP2ckmfA860PUPS8aQprSMBbM8kmQleAny2t7yC/wCul7SGFAb/hNzGjtR2kFoq6U3ADNs9ue4Q4E3A0iI3WfSJuidZw2XWAZMKygZBEAxe2uhSLGk4SaFcbPuyXHwy0Lt/KckLduMuJI/afSVNrCr/LbAX8ELbe9v+Uy5/3PZHa3ThBJLt5VFJ9+XRyaPAcbQ5R/0PgVuzW69JRpwfFJQNgiAYtLTLUK9krJ4OzKuyVy8CXgrcALwMuD/X/yfgb9lQfygwAtgkOYVtA4VCXdmeT7abSNqBlIOlVJisot5fn5d0FfDiXPRO23eUaSgIgmAw0sbV8lOBtwNzJc3JZWeRQs5fkD2z1pAyOAK8AThJ0npgNfDmrEDagu0nACT9wPZJReUKKZW8gnIxKbfx02W26y6ACYIgGPS0cfGj7Zupv6Rlk4yLtr8EfKk9rSckXVFdBBzRG+HY9muaXaPo9Ncv2fDRjQb2Ji3fP7CgfBAEweCkg1fUA0gaDxzNBvfkRcA1tfKskNyH/wx8N9cVMIX6obQ2oZCh3vZzbB+Ut/1IhqKbizYSBEEwGBGdHVBS0knA7aTFktuQ1rMcQVpYWWtKawop/fAngSdt3wCstn2j7RuLtNmSH15e8fm8VmSDIAgGFZ09UvkkaQHkRqOSHBn+T1Q5XGU34q9KujT/fZSSeqKoTaVyZf0Q4FCg/EKCIAiCwYRBPR2tVURttddDg5BkthcAb5L0r6TQLoUpqoHGVex3kWwsM8o0FARBMBjp5Bz1wOeB2/MC9n/ksj2BV5AWSzbE9i9Jz/vCFFUqf7a90Qr6vOqy41bVB0EQ9CsdrFRsfz97dL2SZKgXab3LJ2wXWiFflqIr6j9RsCwIgmCropMN9QBZeVyft98A1/eVQoHmUYpfRYoavJuk/604tS3FMz/2G3uNeYJvPu+HpWT+vm6nUvXvWVU+9tczRjxZWmZNT7l4SQBH7zavtMyqZ5SLS7V4bfn4Yg+sKJSFdCP2HrvJwuC2M3Jox/2En2bdkL6P/TViaHdpmTEjyvfrqeHlYnmtHlE+JtvyoduUlmkbHTxSySlGvkWKbryANFLZXdIy4IO2b293m82mvxYBs4DXkNzMelkBnNnuzgRBEGxRDPAopAAXAe/rjfnVi6TnA9+jRHoRSXNtP6dZvWZRiu8E7pR0se3Ofa0LgiAYAERbk3T1BWOqFQqA7Vty0sWNkHRcnesI2KVIg82mv35m+3jgjlqx920f1KwBSfNJI5tuoMv2lBw6/71scEs+y/avasgeDVwADCVlhPxis/aCIAj6lfaF2+oLrpL0S9J6lF7vrz2Ak4BamXR/ClxM7Um9UUUabDb9dXr+++oiF2vAETUiXX7V9vn1BCQNBb5Ocn1bANwm6Qrbf97MvgRBELSNTp7+sv3v2Tb+WjZ4fy0Avl7rRR64Czjf9t3VJyRV52qpSbPpr4fz7gdtf6yqgS8BH9tUqm0cBvw1pxVG0k9IH0wolSAIOoM2BpTsK2xfBVxVsPoZ1F/s+PoiFyjqUvyKGmWvKihr4FpJsyWdUlF+mqS7JF2YQwZUsxsbhmuQtOsmrleSTpE0S9KsZU+U92YJgiDYHNRTbBuQvknjJX1R0jxJT+RtXi7brrq+7d/Viz5ve1at8moaKhVJH5A0F3hWVgC92wOkYVIRpto+lKSETs1pL78J7AtMBh6mdgTMWiEEatl1ptmeYnvKdjv0vRtmEARBJZ2sVICfkdIAH2F7B9s7kAJKLqOPFq83s6lcQho2/Tfw8YryFbaXFGnA9qL897GcOfKwnPoSAEnfoXau5AUkg1Ivu5NcnIMgCDoD0+mG+kk578rT2H4E+KKkd/ZFgw1HKraftD3f9ltsP0jKLmZgbE7c1RBJYySN690HjgLulvSMimqvBzYxCgG3AftJ2lvSCFJ+5OoEMkEQBANKh6+of1DSRyXt/HR/pZ0lfYyNzQtto5BNRdKxku4HHgBuBOZTzPCzM3CzpDuBW4Ff2r4aOFfSXEl3kYZiZ+Z2dpX0K4C8LuY04BpgHvAz2/eUubkgCII+xwW3geHNwA7AjZKWSFpCiv01ATi+loCk/SW9XNLYqvKjizRYNKDk54DnA7+2fYikI4C3NBPKnlubrNi0/fY69ReRwsL0Hv8KqOX2FgRBMOD0Julqy7WkPUjrSXYhhaafZvuCilAro0jhsT5o+1ZJJ7LBA3cl8IG8YP1pcoyvj1HQU1fSvwOnkl7kp0s63fbl+fQXqL22ZSOKen+tt/0EMETSENvXk4zsQRAEWy928a05XcCHbf8z6SX+VEkHAOcCZ9ueDHwqH0OaOXppXoT+WWBama7Xsam8l5TU63WkbJH/Jal3vWLd/CuVFB2pLMtDoZuAiyU9RgcGlBRmlMoFo3vh6L+Xqn/AyIWl6gOscfngkHPX7NG8UhU/frB8Ms4Va8oF+2uFXbYtleMHgO2HrSots6qnXHDMkUPK/4RHDCmfLHXc8DWlZcqyrqd8v9Z2l5dZ013+t1w+cGX54JDDJqwsLTO/tERt2uXZldcFPpz3V0iax4a88tvmauPJDku2/1AhfgvJmakMZ5Pif1Uy1PbKfP35kg4Hfi5pL9qsVF4LrCHZPk4k3dg5BWWDIAgGLSWmvyZKqlzrMc12zdGFpEnAIaSUv2cA10g6nzS79MIaIu+mhp07261rNkGyeVfziKTJtucA2F4p6dXAhUDTYJJQUKnYfqri8PtFZIIgCAY9BoqnE15se0qzSnlWaAZwhu3lkj4HnGl7hqTjgenAkRX1jyAplRfVuNzOpARd1flTBPxh0+qcRNUsVHaaOknSt5v1HZoHlFxBbb8Fpba8bY1zQRAEWw9t9OySNJykUC62fVkuPpkNcRgvBb5bUf+gfPyqbPeu5kpgbO/Io6qtG6rLcm76mtj+fZF7aBb7a1yj80EQBFs7bfT+EmkUMs/2VypOLQJeSnIFfhlwf66/J3AZ8Hbb99W6pu1312vP9lvb0/ONKW+pC4IgCDbQvhX1U4G3A3Ml9Y4sziJ5ZF0gaRjJtt0bQ/FTpDUo30j6KKUWqbygpLG9hvd6FKlThlAqQRAEreK2en/dTH0Pq+fWqP8e4D1NLnt5VlCXA7N77eOS9iEtPD8e+A7w81b7XU0olSAIghZJix87N/aX7ZdLOgZ4HzA1R4TvAv4C/BI4OccCaxuhVIIgCDaHzk4n3O+RSUKpBEEQbAadPFIZCEKpBEEQtMoWkPmxvwmlEgRB0DJGxRc/bhUUDSgZBEEQ1KJ9ASX7BElvzX9P6I/2BtVIZUnXWH70RK2wOPV51jblHB92HV4d7aA5ZYNcAuwz4rHSMqfuc31pmSXdY5tXquDJrvLB/ta2EOywP1jRNaq0TH8FYRw1tNxvZr8x5X8vrXwvS1v4/llTcg316PJNrFjX94FRa9JGl+I+ZLcc3qVswMmWiJFKEATB5tDBIxVJnyYl5LoEmCDpU33dZiiVIAiCzaGDMz/aPhtYArwNWGK7z6PL97lSkTQ/pw6eUxX2GUkfkWRJE+vIdme5OZIiP30QBB2H7ELbALLQ9k+A8smgWqC/JruPsL24siCnznwF8FADudU521kQBEHnYaB7i/H+6peODuT011eBjxJe3kEQbKGIYqOUAR6pDDpDvYFrJc2WdAqApNeQhmR3NpEdJWmWpFskva5WBUmn5DqzVi/t+7StQRAEGxGG+o3oj+mvqbYXSdoJuE7SvcAngaMKyO6ZZfcBfitpru2/VVbI6TinAex0wA4x6gmCoH/p4DAtts+W9B8kQ/3uts/v6zb7fKRie1H++xgwk5RsZm/gTknzSUOy2yXt0kD276QENYf0dX+DIAgKY1JAySLbwLGoPw31fapUJI2RNK53nzQ6uc32TrYn2Z4ELAAOrQ6/LGl7SSPz/kRSAps/92V/gyAIytLpNhXbF+e/P+6P9vp6+mtnYGbOSjYMuMT21fUqS5oCvD8nn/ln4NuSekjK74u2Q6kEQdBBGHo6f0l9f9KnSiVPWx3cpM6kiv1Z5Exmtv8APKcv+xcEQbBZmI62qQwEnRmUKQiCYEshBiobMaiUyqquEdz++B6lZH6/fu8+6s0Ghg0t/6s7YMKjpWWeM668HW7+mh1K1V+4arvybSzbvrTMkwvHl5YZ9cjQUvXVVboJ1peLvwmA91pdWuZ9B/2uVP2RQ8oHLX1s3balZUYOKf+hlQ2OubJrROk2VqwZoICStC9JV14Q/gNgF5Kqmmb7AkmTgW8Bo0ipgD9o+1ZJ+wPfAw4FPtkfnl1FGFRKJQiCoN9p3/RXF/Bh27dnB6fZkq4DzgXOtn1Vzjd/LnA4KabXvwM11/ANFKFUgiAIWsWG7vbMf9l+GHg476+QNA/YjWS56R1Wjgcql2k8Julf29KBNhFKJQiCYHMoPlKZWBVUd1pevL0JkiaR1uX9CTgDuEbS+SRP2HJJo/qZUCpBEASbQ3Glstj2lGaVJI0FZgBn2F4u6XPAmbZn5Bhe04EjW+5vHxP5VIIgCFrFQI+LbQWQNJykUC62fVkuPhno3b8UOKzdt9FOQqkEQRC0jME9xbYmKK0Snw7Ms/2VilOLSOGtAF4G3N/222gjMf0VBEGwObTP+2sq8HZgrqQ5uews4L3ABZKGAWuA3mjvuwCzSEb8HklnAAfYXt6uDrVCKJUgCIJWMe30/roZUJ3Tz61R/xH6KUdKGUKpBEEQbA4RpmUjQqkEQRC0zMAl4OpUQqkEQRC0iokoxVWEUgmCINgcYqSyEYNKqQwd0sP2o1eVklmxplywwzWrywe7G9JCQMm7up9RWub+ZTuWlnl86bhS9f3IqNJtjH64nu2xPuPLx2AsTdfo8jIjniwvs+7B8g19Y93h5Rsqyb57lA9a+qzxj5WWWb6+3G/moaXlA5A+9UgLkT7bRSiVjRhUSiUIgqBfsXF390D3oqPo88WPkuZLmitpTlXcGyR9RJJzuuBasidLuj9vJ/d1X4MgCErTxhX1g4H+GqkcYXtxZUHOHfAK4KFaApImAJ8GppDMYbMlXWF7aV93NgiCoDAx/bURAxmm5avAR0kKoxavBK6zvSQrkuuAo/urc0EQBE1xzlFfZNtK6A+lYuBaSbMl9YYXeA2w0PadDeR2A/5Rcbwgl22EpFMkzZI0a92yfrDuBkEQVGIX27YS+mP6a6rtRZJ2Aq6TdC/wSeCoJnK1XIY2+WZyPoJpAOP333nr+eaCIOgAwlBfTZ+PVGxXZimbSYq2uTdwp6T5pNg1t+fgaJUsACoTzu9OzngWBEHQEbQ59P1goE+ViqQxOdcyksaQRie32d7J9iTbk0jK49AcHK2Sa4CjJG0vafsse01f9jcIgqA0bQp9P1jo6+mvnYGZKU0Aw4BLbF9dr7KkKcD7bb/H9hJJnwVuy6fPsb2kj/sbBEFQGAPeikYhRehTpWL778DBTepMqtifBbyn4vhC4MK+6l8QBMFmYW9Vo5AixIr6IAiCzSBGKhsjDyJXN0mPAw/WOT0RWFzn3JZG3EtnEvfSmdS7l71slw+YV4Gkq/P1i7DY9qBfazeolEojJM2yPWWg+9EO4l46k7iXzmQw3cuWwECuqA+CIAgGGaFUgiAIgraxNSmVaQPdgTYS99KZxL10JoPpXjqercamEgRBEPQ9W9NIJQiCIOhjQqkEQRAEbWNQKhVJF0p6TNLdFWUTJF2Xs0hel+OJdTx17uVNku6R1JND22wR1LmX8yTdK+kuSTMlbTeQfSxKnXv5bL6POZKulbTrQPaxKLXupeJcw+ysnUad7+Uzkhbm72WOpGMGso+DnUGpVICL2DSh18eB39jeD/hNPt4SuIhN7+Vu4Djgpn7vzeZxEZvey3XAs20fBNwHfKK/O9UiF7HpvZxn+yDbk4ErgU/1e69a4yJqJMBrlp21Q7mI2sn8vmp7ct5+1c992qoYlErF9k1AdfDJ1wLfz/vfB17Xr51qkVr3Ynue7b8MUJdaps69XGu7Kx/eQkpx0PHUuZflFYdjqJ/VtKOo8/8Fmmdn7Tga3EvQTwxKpVKHnW0/DJD/7jTA/Qk25V3AVQPdic1B0ucl/QM4kS1npLIJBbOzbkmclqcmL9xSpr63VLYmpRJ0MJI+CXQBFw90XzYH25+0vQfpPk4b6P60gqRtSNlZt1ilWMU3gX2BycDDwJcHtjuDm61JqTwq6RkA+e9jA9yfICPpZODVwIkePAunLgHeMNCdaJF9KZaddYvA9qO2u233AN8BDhvoPg1mtialcgVwct4/Gbh8APsSZCQdDXwMeI3tVQPdn81B0n4Vh68B7h2ovmwOtucWzM66RdD7Mpl5PcnRJegjBuWKekk/Bg4nhaR+FPg08AvgZ8CeJG+WN20JmSTr3MsS4GvAjsAyYI7tVw5UH4tS514+AYwEnsjVbrH9/gHpYAnq3MsxwLOAHlIKhvfbXjhQfSxKrXuxPb3i/Hxgiu2OD4Vf53s5nDT1ZWA+8L5e+2rQfgalUgmCFLkdTgAABaJJREFUIAgGhq1p+isIgiDoY0KpBEEQBG0jlEoQBEHQNkKpBEEQBG0jlEoQBEHQNkKpBEEQBG0jlMpWhqSVfXDN10j6eN5/naQDWrjGDWXC+Of6f8kxqqrPTaoVxn2wIukdlWH2JV0saYmkNw5kv4Ktk1AqwWZj+wrbX8yHrwNKK5UWOdH2FX3ZgKShfXn9NvEO4GmlYvtEUgSJIOh3QqlspShxnqS7Jc2V9OZcfngeBfw8J8+6WJLyuWNy2c2S/lfSlbn8HZL+T9ILSeFJzsvJkPatHIFImphXZyNptKSf5MixPwVGV/TtKEl/lHS7pEsljS1wP8+VdKekPwKnVpQPzfd5W27rfbl8iKRv5GRnV0r6Ve+bvaT5kj4l6WbgTfk+rpY0W9LvJO2f6+0oaUa+9m2Spubyl1YkhLpD0rgG/f6Pir6dXVH+i9zePZJOqbiXiyq+szNzn6cAF+f2RtdrKwj6g2ED3YFgwDiOFLriYFJIi9sk9Sb9OgQ4EFgE/B6YKmkW8G3gJbYfyOEwNsL2HyRdAVxp++cAWR/V4gPAKtsHSToIuD3Xnwj8J3Ck7ackfQz4EHBOk/v5HvBvtm+UdF5F+buBJ20/T9JI4PeSrgWeC0wCnkNKgzAPuLBCbo3tF+U+/YYUcuV+Sf8CfAN4GXABKfnTzZL2BK4B/hn4CHCq7d9nhbimVoclHQXsRwpwKOAKSS/JOUHeZXtJVhK3SZqR+7ub7Wdn+e1sL5N0GvAR2/+/vXMJjfKK4vjvH3yEGrQg6iIIUkVErBZR2tSgXQRFuhApKuJKSwsKDW505yoLKxQLLaIBFwEJaYmPjYKJ0MbURxRBEx8ogrrooiKYGF8IhuPinqmfw8ykwRljzfnB8H1z7p07916Ge757zp1zLg0zR0FQcUKpjF3qgTYzGyJFcD4NLAUGgYtm9jeApCukxewJcMfM7vrn24Dv3+L7lwO/AJhZn6Q+l39BMp+ddYU0AThfqiFJU4CPzey0iw4Bq/1+JbAw41+YQlrI64F2j1z7j6Q/85r93duuAb4E2jMKcqJfG4D5Gflk35WcBfZKagWO5uayACv9ddnf13jfuoFGSWtdPtPlt4BPJP0KnAA6S81LEIwGoVTGLkW3EMCLzP0Q6XdSqn4pXvLazFqdV1Yo8JyAU2a2cQTfoSJt5cp+MLOON4TS18O0+dSvVcCApwjOpwqoM7PnefIfJZ0gBZjskdRgZoUiFgvYbWbNeX37iqSw6szsmaQuoNrM+iUtAlaRTHzrSYnNguC9IXwqY5duYIPb6aeRdg4XS9S/SXpKnuXvNxSp9xjI+hDukUxNANnTSN2k7IhIWgAsdHkPydw2x8s+kjS31EDMbAB4JKneRZsyxR3AVknjvb25kiYBZ4Bv3LcygxTJtlDbg8BdSev88/KFHdJO4d9EXJI+8+tsDx+/B7gEzCvS9Q5gS85nJKlW0nTSbqrfFco80u4tZxqsMrMjwC5gsbeTP+dBMGqEUhm7HAP6gF7gD2BnqXwZ/jS+DTjpDuz7wKMCVX8DdriDejbwE2lRP0fy3eTYD9S42WsnrtDM7AHpNFObl/VQfFHOshnY54767M7hIHCDlGTqGskvNA44QsoTkpNdKDIeSErqW0m9wHVgjcsbgSXuZL8B5EL2b3dneq/3pWCKZDPrJCXzOi/pKnCYpBxOAuN8/E0+BwC1QJebJFtIaQPw+wPhqA/eByL0ffCfkVRjZk+UnAj7gNtm9vMo9aWLt3ROZ8YzlaTUlv1fE1HlI6mFzIGJIHhXxE4lGAnf+VPydZKJpnmY+pXkIdCiAn9+HAHHfTx/AU0fkEJpBVZQ5NRZEFSS2KkEQYWR9CnpRFqWF2b2+Wj0JwgqSSiVIAiCoGyE+SsIgiAoG6FUgiAIgrIRSiUIgiAoG6FUgiAIgrLxCmFU979ziQn/AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ds.tmpprs.isel(time=0,lev=0).plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Top level = 0.4 mbar, i.e. at the top of the atmosphere."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.QuadMesh at 0x20cc83d8688>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ8AAAEXCAYAAAByAUkhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOydeZxcVbW2nzeBBEggYQogBII4oCAEjDjECQRErqCCIDiAAwYH7gWcx08Q9CKKXPReUSZFBEVFBBEZVAZRGQKEyaCIRAxhMEBCmJP0+/2xd8FJ9amuU9XV3dWV9fA7vzq1zx6rQq3ee6/9LtkmCIIgCIaTMSPdgSAIgmDlI4xPEARBMOyE8QmCIAiGnTA+QRAEwbATxicIgiAYdsL4BEEQBMNOGJ8uRdKmkh6VNHak+xJ0BknTJFnSKiPdlyAYacL4dAmS5knaufbe9t22J9pePpL9aoSk8ZJOlfRPSUsk3SjpTXV53iDpdkmPS7pM0maFZ9+QdEcue7ukA+rKTpd0fS57vaTpTfrTsK38fGdJN0h6TNK/JO3bZGynSXpE0n2SPtYg35eyMdm57PloQNLheYyL85jHVygz6scdjDxhfIJ2WQX4F/A6YBLwReCnkqYBSFoP+EVOXweYDZxdKP8YsEcueyBwgqRX5bLjgPOAHwFrA6cD5+X0fjRrS9KLgbOAz+f2pgPXDzC2I4DnA5sBOwKfkrRbXZtbAG8H7h2gnq5G0huBzwBvAKYBzwWObFJm1I876BJsxzXCF3AG0Ac8ATwKfIr0Y2BglZzncuBo4E85z6+AdYEzgUeA64BphTq3BC4FHgL+Cuw7DOO4Gdg7388C/lR4NiGPb8sGZc8HPp7vdwXuAVR4fjewW4OyA7ZFMjxHtTCOe4BdC++PAn5Sl+c3wO7APGDnivXWf6eTgFNJP+T35O93LDAeWARsXSi7fh7TlA5+X2cBXy28fwNwX5MyLY87rrjKrpj5dAG230P6cd3Daant2AZZ9wPeA2wMbAH8Gfg+6a/9ucCXACRNIBmes4ApwP7AdyRtVVappO9IWtTgurnKGCRtALwAuC0nbQXcVBjjY8CdOb2+7OrAy+rK3my7qP10c1nZim29Irdzi6R7Jf1I0joNxrE28Jxiffl+q0KefYCnbV/YoD9VOR1YBjwP2I5kdA+y/RRpJrd/Ie++wBW2Hyjp86sH+P4WSXp1g/ZX+Nzy/QaS1i3L3MFxB0EYn1HG923faXsx6S/QO23/1vYy4GekHzCANwPzbH/f9jLbNwDnkJZL+mH7I7YnN7i2adYpSauSZmCn2749J08EFtdlXQysWVLFd0k/fBe3UbZK/k1IRntv0nLa6sC3B6irVr5fXZImAl8FDmtQvhLZWL8JOMz2Y9moHE/6AwPSHw5F4/POnNYP21cN8P1Ntn1Vg27Uf261+36fc6fGHQQ1wutmdHF/4f6Jkve1H87NgJdLWlR4vgppea+jSBqT630aOKTw6FFgrbrsawFL6sp/Hdga2LEw02lYVtKmwF9qibYnVmjrCZLh/ltu86vAb/P9d4F353xfBU4slH+ypK4jgTNs38Xg2AxYFbhXUi1tDGkfDeD3wOqSXg7cR9qnOneQbdZT/7nV7peU5O3UuIMAiJlPN9FJefF/kZZoin/9TrT94bLMkr6b3brLrtvKyuRyIu1ZbEDa61laeHwbsG0h7wTSUuFthbQjSX/972r7kbqy26jwqwxsA9zmZ70AJ2bDU6Wtm2nw+dr+UKG+r9p+mLQHs20h27aFut4A/Ff2ELsPmEpytPh0g4+pEf8CngLWK3xHa9neKverD/gpafbzTuAC22VGAUmvGeD7e1TSaxr0YYXPLd/fb/vBkrydGncQJEZ60ymudAFXA7MK76fR3+HgoMLzo4EfFN7vDPw9368J/JO01LRqvl4GvKjDff5u7vfEkmfrk5Zx9gZWA74GXF14/lngDmCjkrLjcv8PJW2+H5Lfj2vQj2ZtvR+4i+TNtQbpR/2MAcZ1DHAFydNuS5Ix2i0/WxfYsHD9C9in9hmQPOUub1Bv/Xd6HnACacYxhmQwX1fI//Lc9q3AW4bg39xupFnVi/NYfw8c0yDvgOOOK65WrxHvQFz5i4C3kJwOFgGfGIzxye9fCPwa+DfwYP5hmd7B/m6W+/ckafmmdr2rrk+3k5a9LmdFbzyT/vIvlv1c4fl2JHfoJ4AbgO2a9KdhW/n5kfmz+DdpmXDtAeoaD5xG8iK8H/jYAHnnUfD6Is0Ev9Igb/13Oom0zDefZDxvBParK/N3ksdiqeHtwPf4sTzGR0jOK+MLz24rfp8DjTuuuFq9ZEcwuSDoFJLmAG9w+dJVEASZMD5BEATBsBMOB0EQBMGwE8YnCIIgGHbC+ARBEATDTk8dMh07YYJXnVyqmtIQD0fAgrGt76uNG7es5TLjxwx9mdXHPN1yG6up9X61w2rqnb+lnnTfMLTR+v/+T/SVarsOyFN9Q/8z004bj/7t/oW21x9Mu2/ccYIffKia8Pz1Nz91se3dmudcOegp47Pq5HXY7EOl6vcNWTpp6P8n75vU+o/vtE3+3XKZzdZ8uOUyW6zRWjtbr/6v5pnq2HJcPzmyIWHLVScMSzvDwe1LHxv6Np6e0nKZW5+Y2nKZOx8f1O97Jf65ZO2Wy1yx83H/HGy7Dz60nGsv3rRS3rEb3bHeYNvrJXrK+ARBEAwnBvoY+j9ge5EwPkEQBG1izNLujPfY9Qz5IrlShM5bJM2RNDunTZd0dS1N0g4Nyh6Yo13eIenAoe5rEARBq/RV/C9YkeGa+exoe2Hh/bHAkbZ/I2n3/P71xQI53sqXgBmk2e31ks53En4MgiAYcYxZHgf122Kk3IPMs/Ltk4AFJXneCFxq+6FscC4lCSEGQRB0DX240hWsyHDMfAxcIsnA92yfRApIdbGkb5AM4KtKym3Ms7FNIIkvblyfSdIsUhhlVpnUusdLEARBuxhYHoalLYbD+My0vUDSFOBSSbeTImoebvscSfuSlIB3riun+oooicmSjdlJAKttPDX+FQRBMKzErKY9hnzZzfaC/PoAKRLjDsCBpBj1kMI/lzkczCcFrKqxCeXLc0EQBCOCgaV2pStYkSE1PpImSFqzdg/sSgqMtQB4Xc62EymoWD0XA7tKWlvS2rnsxUPZ3yAIglYwZnnFK1iRoV522wA4N0dDXgU4y/ZFkh4FTpC0CikY2SwASTOAD9k+yPZDko4Crst1fdn2Q0Pc3yAIguoYlnfIrkiaCvyQFCm2DzjJ9gmSjgA+SAqECCno4oX5iMpJteLAEbbPLan3TJLX8FLgWuBgrxjyfkQYUuNj+x+sGCO+ln4V8NKS9NnAQYX3p5EiSgZBEHQdSeGgYywDPm77hrxidL2kS/Oz421/oy7/rcAM28skbQTcJOlXtuv1vM4E3p3vzyL9xp7YuW63RygcBEEQtI1YXuob1Tq27wXuzfdLJM2lxMO3kP/xwtvVKHHIyvkufKa30rWk/fMRp6eMj5bDuEWtlRm3qLVtr6cnt1Y/AItbVwK+e/FGLZeZN6l1Acdr1tyspfwbTnpBy220I3jaDsMhkjpcbDlu6IVFu5V2RELnzR968dIyDPRVX3Zbr6bykjkpe+v2Q9I0YDvgGmAmcIikA4DZpNnRwznfy0mrQ5sB7ymZ9RTrXBV4D3Bo5R4PIT1lfIIgCIYTA09X99taaHtGs0ySJgLnAIfZfkTSicBRubmjgOOA9wPYvgbYStKLgNMl/cb2kw2q/g5wpe0/VO3wUNI7AVCCIAhGgD6r0lWFPDs5BzjT9i8AbN9ve7ntPuBkSo6m2J4LPAZs3aDeLwHrA63FnBlCYuYTBEHQJknhoDN7PkpuwacCc21/s5C+Ud4PAngbydEASZsD/8oOB5sBLwTmldR7EEmu7A3ZgHUFYXyCIAjaxIjlnVtAmknak7lF0pyc9jlgf0nTSbZuHnBwfvZq4DOSlpKc7j5SE3CWdCFwUD7k/13gn8Cf87GXX9j+cqc63S5hfIIgCAZB1SW1ZuQjKGWVXViShu0zgDMaPNu9cN+Vv/Nd2akgCILRgBFPe+xId2NUEsYnCIKgTdIh0/DbaocwPkEQBIOgUw4HKxthfIIgCNrEFssdM592COMTBEEwCPpi5tMWYXyCIAjaJJ3ziZlPO4TxCYIgaBMjlnanJ3PXE59ai7QqXAptipG2wZjFrX+dT7dYZt6S8S23MY/WRR9XW/Oplstcw9CLpLZDO8Kq57WYv1VR1W7mvsVrjXQXWmJ5h875rGyE8QmCIGiTDiscrFQMufGRNA9YAiwHltmeIelskg4RwGRgke3pVcoOdX+DIAhaoS+83dpiuGY+O9Y0hwBsv6N2L+k4YHHVskEQBN1COBy0z4guu2UV132BnUayH0EQBO1gFHs+bTIcJtvAJZKulzSr7tlrgPtt39FGWQAkzZI0W9Ls5Y+vvNEfgyAYfmxY6lUqXcGKDMcnMtP2AklTgEsl3W77yvxsf+DHbZYFIIehPQlg9Q2nVg9oGwRBMGjUM4dMJW0IYPs+SeuTJgd/tX3bULQ35DOfHE8C2w8A55Kj8ElaBdgLOLvVskEQBN2AgeUeU+nqZiQdDPwZuFrSh4ELgDcDv5D0gaFoc0hnPpImAGNsL8n3uwK1IEY7A7fbnt9G2SAIgq6gRxwODgG2AlYnBZ57Xp4BrQ1cRoqw2lGGetltA+DcHD1vFeAs2xflZ/tRt+Qm6TnAKTkQ0kBlgyAIRhyjjgWTG2GW2n4ceFzSnbbvA7D9sKQh2c4YUuNj+x/Atg2evbckbQGwe7OyQRAE3YChV5wJ+iStansp8B+1REmrMUTbMz0xXwyCIBgZxPKKV9OapKmSLpM0V9Jtkg7N6UdIukfSnHztntN3KKTdJOltDerdXNI1ku6QdLakcSXZ9iLZUuq2QtYFPt7ih1KJnjDZw8lw6bSturj1vwuWTuobgp6sSDv6ce3QquZcO9w9f8KQtwFwNxu1XKZv0rKW8l+zZmu6dgAbTnqk5TLt0KpW25Nt6AcO17/LekxHFQ6WAR+3fYOkNYHrJV2anx1v+xt1+W8FZtheJmkj4CZJv7Jd/4/na7n8TyR9F/gAcOIK47DvLuuQ7XuAewY5rlJi5hMEQTAIOjXzsX2v7Rvy/RJgLrDxAPkfLxia1cgzlyL5IP9OwM9z0unAW1sa4BARxicIgqBNbNHnMZUuYL3agfh8lR6cB5A0DdgOuCYnHSLpZkmnZQ+0Wr6XS7oNuAX4UMmsZ12SdmYtfT4DGLThJIxPEATBIGjhnM9C2zMK10ll9UmaCJwDHGb7EdIS2RbAdOBe4LhaXtvX2N4KeBnw2ewgsEJ1JU10xWH8MD5BEARtkoLJja10VUHSqiTDc6btXwDYvt/2ctt9wMmUHLa3PRd4DNi67tFCYHI+1A+wCbCgpN3dCveTJJ2aZ1pnSdqgUudbJIxPEARBmySHA1W6mpH3Z04F5tr+ZiG96LHyNpKjQc2LbZV8vxkpTM28Ffpnm3RI9O056UDKYxV+tXB/HGmGtQdwHfC9pp1vg/B2C4IgGAQdVDiYCbwHuEXSnJz2OWB/SdNJtm4ecHB+9mrgM5KWAn3AR2rhZyRdCByUz05+GviJpKOBG2muVjCjEF/teEkHdmR0dYTxCYIgaJNOKhzYvoryPZoLG+Q/AzijwbPdC/f/oLku5hRJH8vtryVJedYEQ7RCFsYnCIJgEPT1xu7FycCa+f50YD3g31npek7DUoMgjE8QBEGb2PREMDnbRzZIvw84YCjaDOMTBEHQJkYs66vmyRasSBifIAiCQVBFvSDoz4DGR1Iz8ScB99p+Qee6FARBMDqouVoHrdNs5nOn7e0GyiDpxg72Z1CMWQbjF7V2ePepyb3zD6cdMdKVmXGLRroHA1EmPNyYpye3lh/g3mf2l6szLOK1bZQZuX/76qSwaFciacNafJ/8fvsKxZbavmWgDM2Mz94VGqmSJwiCoCfp6/1lt1MpxPgBriAdPh1o4JsD0waqdEDjk/3DB6RKniAIgl7EhqU97nBg+z/qkq6zvdNAZST9vlm9leaLkvbKgYgWS3pE0pIK+0G1svMk3ZIDHs3OaWcXgiDNK5zmrS+7m6S/Svq7pM9UaS8IgmC4qB0y7YS8zmihmeGpmqfqYuWxwJ62J9ley/aatluJELWj7em2Z+SOvSO/n04S0ftFfQFJY4H/A94EvJgkMfHiFtoMgiAYcvpQpWu0IumCBulfrns/VtKZVeutanzuz6qpHSUL6e0L/Ljk8Q7A323/w/bTwE+At3S6D0EQBO3SSWHRLuaDDdI3lfRZAEnjgXOBO6pW2szVeq98O1vS2cAvgadqz2uS300wcIkkA9+ri2HxGpJhK+vwxsC/Cu/nAy8v6eMsYBbAuAlr1z8OgiAYUnrN203SOiRB7IdJN/c2yPo+4MxsgHYEfmP7+KrtNPN226Nw/ziwa+G9KVkuK2Gm7QWSpgCXSrrd9pX52f6Uz3qgYhCkbMxOApiw3tSuCJIUBMFKwuif1QAgaVPS9sobgEUpSWsBvwc+Y3teIW/R1foEUsiFPwJXSNq+Fgq8Gc283d6XG5tp+491nZ1ZpYEs6Y3tBySdS1pOuzLHodgLeGmDovOBqYX3pUGQgiAIRgoDy3pj5nM28D/Au2wvh2f23fchbXm8opD3uLqyD5P25Y8jfSRNnQ2gurzOt4H6g0VlaSsgaQIwxvaSfL8rUNuk2hm43fb8BsWvA54vaXPgHmA/4J0V+xsEQTDk9JDCwXq2zy4mZCP0E0lH1aXv2IkGm+35vBJ4FbB+jvVQYy2ginP7BsC5ya+AVYCzbF+Un+1H3ZKbpOcAp9je3fYySYcAF+e2TrN9W4U2gyAIho0eMT7XS/oOKZxCba99Kiny6QoqNpLebLvUA66VPM1mPuOAiTlfUYvjEZ4Ny9qQfAB12wbP3luStgAoBkG6kAaBlIIgCEaaTgaTG2EOAD4AHEly9hJp6+N8+kc+/bqkexhY4eCrQPvGx/YVkq4CXtIo3kMQBMHKzGg+w1MjH2c5MV/NuB/4ZpM8TV2um+752F6eXe+6Hi2H1Ra1KnzY2mbh+DbEKHtJvLSbaVVUtpdo599lewz9v+VR9f+LO7fsJmkq8ENgQ6APOMn2CZKOIJ21+XfO+jnbF0raBTiGtEL1NPBJ2/1kbSRNB74LrAYsAz5i+9q6PKuQZj5vJc18THLwOg841fbSZ4Zsv74T463qcHCjpPOBnwGPFTpRxdU6CIKgJzGwrK9j3m7LgI/bvkHSmqR9mEvzs+Ntf6Mu/0Jgj3yUZWvS/vjGJfUeCxxp+zeSds/vX1+X5wySi/WRpOU2SB7GBwI/At4xuKH1p6rxWQd4kBVd6Kqe8wmCIOhJOrnnkw9z3pvvl0iaS7kxqeUvOgLcBqwmabztp+qzkpzEACZRfmRle9svrEubD1wt6W8tDKMylYxP7bxPEARBsCKubnzWq4krZ06qU3x5BknTgO2Aa4CZwCGSDgBmk2ZHD9cV2Ru4scTwABwGXCzpG6R9hleV5HlY0j7AObb7ch/GkM751LfVEaqqWm8i6VxJD0i6X9I5kjYZig4FQRCMJloQFl1oe0bhamR4JpIElw+z/QjJCWALYDppZnRcXf6tgK8BBzfo4oeBw21PBQ6nv/capKMvbwful/Q3SXcA95GEAPZr0M81JH1R0sn5/fMlvbnhB1VH1cXK75Nc7p5Dmgb+KqcFQRCstNidFRaVtCrJ8JxZ21O3fb/t5XlGcjJJJaaWfxOSoOcBtu9sUO2BPLtF8rNi+WfH4Xk52sD6wCuBV9qektPualDv90lan6/M7+cDR1caKNWNz/q2v297Wb5+AKxftZEgCILeRCzvG1PpalpTOo1/KjDX9jcL6RsVsr0NuDWnTwZ+DXy2Xv6sjgXA6/L9TjRwg5Y0UdLbSWd+9s/x1Abq+Ba2jwWWAth+ghbcIasan4WS3p3jNYyV9G6SA0IQBMFKja1KVwVmAu8BdioE29wdODYH5LyZpB59eM5/CPA84IuF/FMAJJ0iaUbO90HgOEk3kQ5/zqpvWNK+wGXAbrneHXJf5kh6SYP+Pi1pdbLgs6QtKEQ9aEZVb7f3A/8LHJ8b+lNOC4IgWGnppLab7asonzmUqrzYPpoGy1y2D6qrt5GAc40vAK+w/bik9UjLfm+UtA1JtbrMSeFLwEXA1BxEbibw3ibtPENVb7e7gT2rVhoEQbBS4LTv0wMIeCLfPwZMAbB9cw6tsGLmtER4O8kh4RW5/KG2F1ZtsJLxkbQ+aeo2rVjGdsx+giBYqekFeR3S7OoiSVcAbyI5JtQCy/UboG1L+qXtl5L2nVqm6rLbecAfgN8Cy9tpKAiCoNcwLZ3z6VpsfzrvL70Y+LLtmrLCIhqHzrla0stsX9dOm1WNzxq2P91OA0EQBL2LWN43+o0PlEcRyO7djZwIdgQOlvRP0lKdUhFvU6W9qsbnAkm7584FLbIyC15CO2Kvw8P4h5Y2zzRKeGqdVUe6Cx1kdEUG7YWZT5u8aTCFqxqfQ4HPSXqK5NNds3D9NqKCIAhWFuyV2vgM6q/qqt5uaw70XNJWjaKMSpoHLCHtFS2zPSOn/yfJn3wZ8Gvbn6paNgiCoFvokWBy7fBrkgESKVzD5sBfga2qFK4682nGGTTelALYseiCJ2lH4C3ANrafqh2MqlI2CIKgm+gFV2tJu9m+KN9PIgWLexlJTeFw2/fXl7H9kro6tqexvlw/OrW42qrp/zBwTE2B1fYDHepHEATBsGFEX9+YSleX89XC/XEkAdM9gOtIh0ybYvsGksGqRKdmPgPZfgOXSDLwvazk+gLgNZK+AjwJfKKBu15Z2SAIgq6hByY+9cywPT3fHy/pwLJMkj5WeDuGtPr177K8ZXTK+AzEzBxpbwpwqaTbc7trk07Gvgz4qaTn2v0msP3K2r6ymEHSLLJW0fjVJw/5YIIgCJ6hdxwOpmRjImAtSSr8HjeathV9AZaR9oDOqdpgp4zP040e2F6QXx+QdC5JsG4+8Is8uGsl9QHrUWc1G5S9si7PScBJABPXntqDf4QEQdDV9Mavzsk8a0xOJ/8eS9oQmNOgzF9s/6yYkAPS/axB/hWoKq8zE5hj+7GsaL09cILtfwLYfkWDchOAMTkk7ARgV+DLwKMkae/LJb0AGEeKR16lbBAEQdfQCzMf20c2SL+PFGKhjM/S39CUpZVSdeZzIrCtpG2BT5FiTvyQZ2NENGID4NykQccqwFm2L5I0DjhN0q2kWdOBWSvoOcAptndvVLZif4MgCIaFXvB2A5C0JSlY6DW2Hy2k71b87ZX0JmB3YGNJ3ypUsRZp+a0SVY3Psmwc3kKa8ZzaaBOqiO1/ANuWpD8NvLskfQFpUA3LBkEQdAs2uPs92Zoi6b+AjwJzgVMlHWr7vPz4q6TQCTUWALNJkQ6uL6Qv4dlYQ02panyWSPosKbjQaySNBXpJzyMIgqAtemTm80HgpbYflTQN+LmkabZPoO4oje2bgJsknWW7bY2qqsbnHcA7gffbvk/SpsDX2210NPPU5OFZ313Z9eDa0V1b9cHHhqAnK6J/PzzkbQB4/bVbyr9qG3GFl647oeUyvaUh1yF643/VsbWlNtvzJL2eZIA2o/E5zmmS/pukhL1aLdH2c6s0WGm+mDedzgHG56SFwLlVygZBEPQu1UJojwKnhPsk1c72kA3Rm0leb43CaH+f5A+wjKRw/UOS2k0lKhkfSR8Efs6zJ103Bn5ZtZEgCIKexRWvJkiaKukySXMl3Sbp0Jx+hKR7JM3J1+45fRdJ10u6Jb/uNEDd/ynpr7neY0uyHADct8Kw7GW2DwBe26Da1W3/DpDtf9o+guTFXImqy24fJZ2xuSZ36o4memxBEAS9T2cPmS4DPm77BklrAtdLqgV1O972N+ryLwT2yAfxtwYuJk0MVqCKlqbt+Tnv+sAmuS932X7U9h8b9PdJSWOAOyQdAtxDDr9dharG5ynbT2e3ZyStQq+sdAZBEAyGDhkf2/eSNNXI5xvnUmJMCvlvLLy9DVhN0viaZmaBplqakl4MfAuYBmwK3EhSPbgCONT24pIuHAasAfwXcBRp6a2pF3SNqj6CV0j6HLC6pF1Ih4h+VbWRIAiCnqVDy25FssfZduTVJuAQSTdLOk1SmTfK3sCNJYYHntXSvEbSFZLKxD9PAz5q+3nAq4HbbW8O/JF0rrO+f2OBffPMaL7t99ne2/bVVcdY1fh8hiR9cwtJMvtC4AtVGwmCIOhZqhuf9STNLlyzyqqTNJHk4HWY7UdIm/pbANNJM6Pj6vJvBXyNxuEMilqanyRpadZP11a3/VcA29eSnQxsn0zyZltxyPZy4KUl9VSmajC5Pkk/Aq6sdTAIgmClx7Sy7LawWUBMSauSDM+Ztn8BUIylI+lk4ILC+01InscH2L6zQbVVtDTvlPRF4HfAXmQ9t9yfRnbiRuA8ST8DnjnnUOt3M6p6u+2ZO1MLNjRd0vlVygZBEPQyKZR286sZeRZxKjDX9jcL6RsVsr2NFOANSZNJStKfHcApAJJn8k65TKmWJvB+krDo54CngENz+ho01nZbB3gw171Hvt484CALVHU4+BLJ2+1yANtz8ppkEATByk1fx7zdZpJUZG6RVFOS/hywfz6DY2Aezy6vHQI8D/hinrUA7JqjAJwCfNf2bNJ+Tj8tzWLDtheRdDupS18MlO7j2H5f2yOlNW23xYNY3guCIOhJ1CG/X9tXUa4mcGGD/EcDRzd4dlDhvlRLc7DkWdSJwAa2t5a0DbBn7ldTqjoc3CrpncBYSc+X9G3gT+11OQiCoEeo6mzQmwdTTiaFUFgKYPtmYL+qhasan/8EtiKtBZ4FLCb5eAdBEKzEKDkcVLl6jzWyZ1yRzoVUyP7cR9r+JPD5Fjs3rHgsPDl59MubQ28JmLb3nbQjYNmaUOZwCJECLLu/35m+pgxHfPt2PuGn1mk9VH2r3387//afbr1bnaMHZjVZOOADJIeG55BGtQA4Dzi1gXr1Qklb5LxIejv5kGwVmv4bt71c0kurVhgEQbBS0QPGhyQIugg4guSaDUlm560FY7sAACAASURBVEDgR6TIBvV8FDgJ2FLSPcBdwLuqNlj1D6wbs2t1W/7cQRAEPYnppLfbSLK97RfWpc0Hrpb0t7ICOeDnzpImAGNsL2mlwarz4bb9uSXNy6qrcyTNLqQ3U1lF0m45z98lfaZiX4MgCIYNudrV5TwsaZ8sFAqApDGS3gGUBrGStG4Oo/0H4HJJJ0hat2qDVRUOBuXPDexo+5lDTVVUVvNe0/8Bu5As8HWSzrf9l0H2JQiCoHN0v2Gpwn4kiZ7vSHqY5PI9CbiMxh5sPwGuJOnKQVpyOxvYuUqDlYxPtm71LAZmF+J8t0JTlVXSoda/56kdkn5CMlhhfIIgCDqI7XnkfZ08e1FxwtCAdWwfVXh/tKS3Vm2z6rLbaiRRuzvytQ1pKe4Dkv6nSVkDl+RgRzUhvSoqqxsD/yq8n095rIpZNaG+ZU8Oj/dSEARBjR5ZdnsG2w8Ca0raS9KWA2S9TNJ+eXlujKR9SXI/lajqcPA8YCfbywAknQhcQloSu6VJ2Zk52NEU4FJJt7OiyurLSCqrz62TfCjbxev3Fdo+ieRxwYT1po6irzgIgp6gB87wSPql7bfm+7cA/0OSUztG0ldt/6Ck2MHAx0jecJAmM49J+hhg22sN1GZV47Mx6RBFLaDQBOA52Q27LH7EM9hekF8fkHQuaTmtisrqfGBq4f0mJL/zIAiC7sBA30h3oiNsVrj/NGmycZek9UhK1z+oL2B7zcE0WNX4HAvMkXQ5aUbyWuCr2cXut40KFV3w8v2uwJeBR0mec5cPoLJ6HfB8SZuTwrPuB7yz6sCCIAiGg9G0pDYAxVGsYvsuANsL8+SglKznNo2CLal6BKeqt9upki4kzVoEfK42oyEFJ2rEBsC5WZB0FeAs2xdJGkeJyqqk5wCn2N7d9jKluOAXA2OB02zfVqW/QRAEw0ZvGJ9tJT1C+n0fL2lD2/fl3+qxZQUknUba/7+NZ+d/BjpnfHKciTcAz7X9ZUmbStqhRNdnBbKn2rYl6aUqq9mg7V54fyENFF2DIAi6gh4wPrZLDQwpnk+jCKmvsN0vymlVqnq7fQd4JbB/fr+EdAYnCIJgpaWqp9toWpqTtLakNSHF+bH95wZZ/yypbeNTdc/n5ba3l3Rj7tDDeTrWVfSt0rooYauChOMWtZa/2xkOAdN2xEuHQ4y0HZHM8eu2Jl4KsOr6a7dcplt/q8Y/VKYvOTBPTh4/BD3pInpAXidveRxDOks5Ebgnb5ecBnylgbDo6SQDdB8p4oFIXm7bVGmz6v/hS7PiQE29dH16xccjCIJgEPTIzOdHpH31ScA+wDnAi0gTlEarXKeRIq/uxrOSa3tUbbDqzOdbwLnAFElfAd4OfKFqI0EQBD1L9xuWKqxr+3JI3mqSPm/7MeAL+WxmGXfbPr/dBqt6u50p6XqS04GAt9qe226jQRAEPcHomNVU4d+S3g38nqTVNg+ecTZrtEJ2u6SzgF+Rlt2A6q7WAy67SVqndgEPAD8mRTK9P6cFQRCs3HQojLakqZIukzQ3q/0fmtOPkHRPjgwwR9LuOX2XLFt2S37dqUn9n5DkfHC0nvcDe5KUa14OHJLT1yGFyi5jdZLR2ZUWox1A85nP9aSPTcCmJGltAZOBu4HNqzYUBEHQk3Ru5rMM+LjtG7K32fWSLs3Pjrf9jbr8C4E9snzZ1qQzkf30LyEZNpIc2t1lz23fDexbkv4gaf+nrMygoh0MOPOxvbnt55IGtYft9WyvS7JuEUguCIKVnk45HNi+1/YN+X4JMJcGxiTnubFw2P82YDVJjVwLjwc+RQdNpaQXSPpdFgtA0jaSKvsCVPV2e1k+8AmA7d8Ar2utq0EQBD1I9WW39WoK/PmaVV4hSJoGbAdck5MOkXSzpNMklfnu7w3cWAtTU1fXnsA9tm9qb4ANOZm0JLcUwPbNNI7904+q3m4Ls0X7EeljfDcpsmkQBMHKS2sOBwttz2iWSdJE0lLXYbYfyVEEjkqtcRRwHGmPppZ/K1IguF1L6loD+HzZsw6whu1r83mgGsuqFq4689kfWJ/kbn1uvt9/wBJBEAQrAx1yOACQtCrJ8JxZ8xqzfb/t5bb7SLONHQr5NyH9Jh9g+86SKrcg7c3fJGkeKTrADZI2bNKPV0v6mKSBjNZCSVvURifp7cC91UZa3dX6IeDQqpUGQRCsNHRoFyW7NZ8KzLX9zUL6RrZrP+pvA2p7LJNJwds+a/uPpV2zbwGmFOqaB8yoj1Iq6VrbO+T7DwIfJRm1L0na3vYxJdV/lBRLbUtJ9wB3kUJpV6KZq/URzSqokicIgqAXER1VOJhJUgzYqc6t+tjsTn0zsCNweM5/CCnQ5xcL+acASDpFUtMlvgJFbapZwC62jyQt1zUyKLa9M2klbEvbr6b6alrTmc9BWWa7ESJtMB1RtcEgCIKewdA42k2LVdlXUR7BuVTZ3/bRwNENnh3UIH1ag+bHZEeGMYBs/zvnf0xSo32cc4DtsxJCjZ8DL22QfwWaGZ+TgWbR6k6u0tBw0DcOHttsaCXnlk4a0uqfYdXF7QhrDj3tCKsOh3hpYug/s3ZEMldbpzXB026mPcHXoWfppBGUmuwNhYNJpHOdAlyI5zOROoMoaUtgK2CSpL0Kj9YCVqva4IDGJ0+7giAIgkb0gPEZYEbUR9pnKvJC0lnPyawoJLoE+GDVNqu6WgdBEAQl9IK2m6TJtvuta9h+nORIUEw7DzhP0isHiPXTlCGfQ0ualzfL5kiandNKtYqqlA2CIOgqOuhqPYIslPRbSR/IXnRNGYzhgeGb+exY79pHuVZR1bJBEAQjTwcdDkaYucD/kM5vHivpKpKQ9Hm2nxiKBivNfAar4RMEQdCz9MbMZ6ntC2y/i3QQ9UyS0Oj8HDah41RddhuMho+BS7Lkd1HLqJlW0UBln0HSrJpW0vJHHyvLEgRBMGT0SCTTZzzabD9h+6e29wJqwtL9C0gzJB0u6euSvixp31ZC7VQ1PmvYvrYuraqGz0zb2wNvAj4q6bXAiSTZh+kkOYbjWii7ArZPsj3D9oyxEydU7FIQBEGH6I2Zz5llibYX2z69mCbpvZJuIE1IVgf+Sor39mrgUkmnS9q0WYOtCIu2peFTk/y2/YCkc4EdbF9ZGMjJwAVVywJXluUNgiAYdkaHYWlKxf33GhNIE4PSvSBJ04Hn0yB2UI2qxqdMw+fdzQpJmgCMsb0k3+8KfLmRVlGVshX7GwRBMOSIckmC0YakMcB7SaEZNiGtbN0BfNf25cW8tv9voLpsz6nSZlVh0X8AOxcNQpVywAbAuVlyexXgLNsXSTojW0eTYoUfDCDpOcAptndvVLZiu0EQBMNCj3i7nQr8E/hv4O3AI8AfgC9Ieontb9cySloF+ABp4vAc0u/4AuA84FTbS6s0OKDxkfSxBukAFJVXy8hGa9uS9Pc0yL8A2H2gskEQBF1FDyy7AS8thMW+StLVtv+fpCuBOcC3C3nPABaRND3n57RNgANJMd/eUaXBZjOfmq7bC4GXAefn93sQey9BEAS9YnyWStrC9p2StgeeBrD9lNTPV2972y+sS5sPXC3pb1UbrKTtJumS3OCS/P4I4GdVGxk2xpq+SZUD6XU1T7UhYDpm8dCfGW5HWHX4RFK7c/V9ybSxI92FUcWIioS2yuhwo67CJ4HLJD1JCq+wH4Ck9envEPawpH2Ac3KAu9qe0T7Aw1UbrPprtSnZEmaeBqZVbSQIgqBn6QHjY/v3kjYD1i0qyuTQCp+qy74fKWz3dyTVjM3awO+pfv6zsvE5A7g2uzubtNH0w6qNBEEQ9Co94nCAbQNNpcxszyPv60halxT/p2UJtKrebl+R9BvgNTnpfbZvbLWxIAiCXqNHlt3awvaDAJJ+aPuAVspWMj75tOpCUkzvZ9JsD3iIKAiCoKfpkUOmrSDp/PokYMeaGrbtPavUU3XZ7dc8+xGvDmxOklTYqmL5IAiC3qRDxkfSVNJ2xoakIG4n2T4hO3h9EPh3zvo52xdK2gU4BhhH2of/pO3fl9T7dZKH8tPAnaSVq36xeyRNAnYDNubZszsXl+TdBPgLcErOJ2AGjWXSSqnkhmT7Jba3ydfzSTI3V7XSUBAEQa8hOiosugz4uO0XAa8g6Vm+OD873vb0fF2Y0xYCe9h+CemMzRkN6r0U2Nr2NsDfSJpsK45DOgC4AXg9sAZJQmdH4Pr8rMgMUsjtzwOLswLCE7avsH1FpZHSZjwf2zdIelk7ZYMgCHqKDs18suTYvfl+iaS5pFlIo/zFfffbgNUkjbf9VF2+SwpvryYpGNTzedJB0xVmOTniwDUUHMyye/Xxkn6WX++nDVtSdc+nqHQwBtieZ6eAQRAEKycG9VW2PuvVRWQ+yfZJZRklTQO2I/3wzySFoDkAmE2aHdWfp9kbuLHe8JTwfuDssiYpN6N9NDhAZ3s+sI+k/yDJ8bREVWu1ZuF+GWkP6JxWGwuCIOg1WvB2W2h7RtP6pImk39fDbD8i6UTgKJJxOIq0t/L+Qv6tSOdudm1S7+dJv99l4RO+AtyQBQX+ldM2BXbJbTbE9q9JNqElqhqfv9heQdEgn3DtPpWDIAiC4aSD3m6SViUZnjNt/wLA9v2F5yuEoJG0CckL+QDbdw5Q74HAm4E35PM8K2D79OzF9kbSUp+Ay4HPlsyyOkJV3ZN+G1QN0oIgCFYqOuVwoKTYfCowtyjaLGmjQrZnQtBk1+ZfkwzEHweodzfg08Ceth9vlC8bmcvy9TvgsqEyPNBc1fpNJJXpjSV9q/BoLapHMu1qVluz2RLpyPDkkvEtl+lWXbul7fm1tN5OG7pzvUK3fvfQ+v9j44aoH0NG52Y+M4H3ALdIqsXE+Rywf1kIGuAQ4HnAFyV9MaftmoNvnkKKxTMb+F9gPCnKKMDVtj9UbDjX/11gEkkkVMAmkhYBH7F9Q8dGmWn2q7CAtMG1J8m1rsYS4PBOdyYIgmBU0UFhUdtXUb65f2FJGraPBo5u8Oygwv3zKjT/A+Bg29cUEyW9Avg+FcPbSLolu343pZmq9U3ATZLOtN29f1oFQRCMAKJntN0m1BseANtX5yCizyBprwZ1iHRAthLNlt1+antf4MaSmA7kQ0sDImkeaaa0HFhme0ajE7slZXcDTgDGkiKcHtOsvSAIgmGl//79aOQ3kn5NOs9T83abChwA1EeQPpvkMVc28NWqNths2e3Q/PrmqhU2YMcS1dPjbX+jUQFJY4H/I7n6zQeuk3S+7b8Msi9BEAQdoxeERW3/V97jfwvPervNB/6vZGJwM/AN27fW1yNp56ptNlt2uzfffsT2p+sa+RrJg2Ko2AH4ew6njaSfkD6YMD5BEHQHPSQsavs3wG8qZD2MxodK31a1vaqu1ruUpL2pYlkDl0i6XtKsQvohkm6WdFqWcKhnY56d/kGywv2kJiTNkjRb0uzlSx6r2KUgCILOoL5qVzcjaZKkYyTNlfRgvubmtMnFvLb/0CiiQfauq8SAxkfShyXdArwwG4radRdp6lWFmba3Jxmrj0p6LXAisAUwnaRlVKaGWub1UbbvdJLtGbZnjF1zQkmRIAiCoaMXjA/wU1II7B1tr2t7XZKw6CKGSEyg2Z7PWaRp2H8DnymkL7H9UJUGbC/Irw/kSKg72L6y9rz+xG6B+aQNrxqbkFy/gyAIugPTKw4H02x/rZhg+z7gGEnvG4oGB5z52F5se57t/W3/E3iC9HFPzAHmBkTSBElr1u5J2kO3NjqxW8d1wPMlbS5pHCk2eH0QoyAIghGlgyEVRpJ/SvqUpA1qCZI2kPRpVtz+6BiV9nwk7SHpDuAu4ArSKdsqG1MbAFdJugm4Fvi17YuAYyXdIulm0tTu8NzOcyRdCJDPFR0CXAzMBX5q+7ZWBhcEQTDkuOLV3bwDWBe4QtJDkh4iabutA+xbn1nSlpLekEVQi+m7VW2wqu7J0aTgRr+1vZ2kHYH9mxXKnmr9Tsbafk+D/AtIcj619xfS4HRvEATBSFMLJjfayRpun6aCB7Ok/wI+SpoUnCrpUNvn5cdfpf+5oFKqersttf0gMEbSGNuXkZwFgiAIVl7s6tcopWTP54OkwHNvJUU+/aKk2pnQ0tg/ZVSd+SzK06srgTMlPUA3CosuF2MWtyZi+eQQdWWwtDqOdhkOQcpuFr1slW4Vom2HDSe1HP+rLTZbc8iEkQfFXztUzyjwZBssR5L03WqMtf0ogO15kl4P/FzSZgyB8XkL6Xf6cOBdJOXTL1dtJAiCoFfphWW3vP9e+oi0d1/kPknTbc8BsP2opDcDpwGVREWhovGxXTy9eXrVyoMgCHoaA9XDaHczG5ACydVPUwX8qS7tAOpWvrKD2AGSvle1wWbCokso99NQas9rVW0oCIKgJ+kJ28MFwMTabKaIpMuL723Pb1TJQEHt6mmm7bZm1YqCIAhWRnph2c32BwZ49s6haLOqt1sQBEFQRg94u9Wf12k3TyuE8QmCIGgX94y223mSjpP02mLwOEnPlfQBSRcDlQ+QViGMTxAEQZukQ6audDWtS5oq6bKsJn1b7eyMpCMk3SNpTr52z+m75GgBt+TXnRrUu46kSyXdkV/7RRGw/Qbgd8DBwG2SFkt6EPgRKTrpgbZ/3u7nVMbwHCYJgiDoVTo3q1kGfNz2DVkT83pJl+ZnZcE3FwJ72F4gaWuSFFm/sDMkUejf2T5G0mfy+35KBsOtKBMznyAIgkHQqZmP7Xtt35Dvl5Dka8qMSS3/jbWoAcBtwGqSxpdkfQvPHpE5HXhrC8MbMsL4BEEQtEtVUdFke9arBb7M16zSOgFJ04DtgGtyUrPgm3sDN9ouk+DYoBaVOr9OaXWYQ0EsuwVBELSNUfVDpgttz2iWKXuVnQMcZvsRSScCR5FM2FGk4JvvL+TfCvgaKWTNqCFmPkEQBIOhg67WklYlGZ4zbf8iVe/7bS+33QecDOxQyL8JcC5wgO07G1R7fy2GWn59oEHb78yv+1Xq7CBZ6Wc+rQp4rrq4x+z14nFD3sTSScPjZ9qqgGk7IqEruxjnFmv8e8jb2Hr1IYld1o8fdqISd86NWpKAU4G5tr9ZSN+otmxGIfimpMnAr4HPNlEWOB84EDgmv57XIN/GkvYlRY0ecnrslzQIgmCY6dzMZybwHmCnOrfq0uCbpGCbzyOFNKjlnwIg6RRJtSW+Y4BdckDQXfL7FZD0JVLguLOAdST9v7Y/j4qs9DOfIAiCQdEh8QLbV1EekqDU/dn20aRAn2XPDircPwi8oUnbR0r6JPBuYJMSt+6OM+QzH0nzstWeI2l23bNPSLKk9RqUXV6w6OcPdV+DIAhapVOu1l3APbZ/AtwzHI0N18xnR9sLiwmSppKmgHcPUO4J2xExNQiC7sTA8lFhWFphWAY0kns+xwOfolcEyYMgWOkQ1WY9o2Tm03MOBwYuydpDswAk7Uma4t3UpOxq+TDW1ZJKT+VKmlU7tLX80cfKsgRBEAwdvaFq3ZMOBzOz9tAU4FJJtwOfp9qBqE1z2ecCv5d0S70vu+2TgJMAxm86tbu/4SAIeo8uNyxV6EmHg5r2kO0HSIehXgdsDtwkaR5pineDpA0HKPsP4HKS3EQQBEF3YJKwaJWr+1kwnA4HQ2p8JE3I6qzkGBG7AtfZnmJ7mu1pwHxge9v31ZVduyaSl73hZgJ/Gcr+BkEQtEqv7PnYPjO//ng42hvqZbcNgHPTwV1WAc6yfVGjzPlQ1Ieyj/qLgO9J6iMZyWNsh/EJgqCLMPSNjmlNtzGkxicvl23bJM+0wv1s4KB8/yfgJUPZvyAIgkFhemLPZyQIhYMgCILBEBOftugp46PlvSP8OW7RSPegk7T+nQyXGOlw0I5I6HAIeA4XwyEUuuW4UqHmYWE07Od0Iz1lfIIgCIadMD5tEcYnCIKgXWxY3juz9OEkjE8QBMFgiJlPW4TxCYIgGAxhfNoijE8QBEG7GOgL49MOYXyCIAjaxuDY82mHMD5BEASDIZbd2qI3DsUEQRCMBCZ5u1W5miBpqqTLJM2VdJukQ3P6EZLuKUR13j2nr5vzPyrpfweod3oOSzMnh5/ZoVPDHwwx8wmCIBgMnZv5LAM+bvuGLMh8vaRL87PjS8IcPAl8Edg6X404FjjS9m+y4ToWeH2nOt0uYXyCIAjapnOB4mzfC9yb75dImgtsPED+x4CrJD2veSdZK99PAhZ0oLuDJoxPEARBu5hWVK3XkzS78P6kHAyzH5KmkeKXXUMKJ3OIpAOA2aTZUSuaTYcBF0v6Bmmr5VUtlB0yYs8nCIJgMFQPo73Q9ozC1cjwTATOAQ6z/QhwIrAFMJ00MzquxR5+GDjc9lTgcODUNkfaUXpq5jP2KVhzXvd5njw1WSPdhYY8Pbm1/MMl+Nk3admwtNMrtCPeOVxinLc/PaVL2+jQ6lMHvd0krUoyPGfa/kWq3vcXnp8MXNBitQcCh+b7nwGndKCrgyZmPkEQBO1i4+XLK13NUIq6eSow1/Y3C+kbFbK9Dbi1xV4uAF6X73cC7mix/JAw5DMfSfOAJcByYJntGYVnnwC+Dqxve2FJ2QOBL+S3R9s+faj7GwRB0BKdUziYCbwHuEXSnJz2OWB/SdNJO0zzgINrBfLv61rAOElvBXa1/RdJpwDfzQE6PwicIGkVkofcrE51eDAM17LbjvXGRdJUYBfg7rICktYBvgTMIH3o10s6v8WNtiAIgqGlc95uVwFla/QXDlBmWoP0g+rqfelg+9dpRnLZ7XjgUyTDUsYbgUttP5QNzqXAbsPVuSAIgqbYydutyhWswHAYHwOXSLpe0iwASXsC99i+aYByGwPFXdT5lPi8S5qVT+3OXvbkY53sdxAEQXOqe7sFBYZj2W2m7QWSpgCXSrod+Dywa5NyZdPPft9gdlc8CWDCelPjGw6CYBhxJWeCoD9DPvOxvSC/PgCcS/K62By4KW+WbQLcIGnDuqLzgamF95vQJSdzgyAIgGdDKlS5ghUYUuMjaULWKELSBNJs5zrbU2xPy5tl84Htbd9XV/xiYFdJa0taO5e9eCj7GwRB0DLuq3YFKzDUy24bAOcm93VWAc6yfVGjzJJmAB+yfZDthyQdBVyXH3/Z9kND3N8gCILKGHDMatpiSI2P7X8A2zbJM61wPxsougieBpw2VP0LgiAYFI5gcu3SU/I6QRAEw03MfNpD7iEXQEn/Bv7Z4PF6QD8VhVFKjKU7ibF0J43Gspnt9QdTsaSLcv1VWGg7zipmesr4DISk2UVpn9FMjKU7ibF0J700ll4ihEWDIAiCYSeMTxAEQTDsrEzGpzRw0yglxtKdxFi6k14aS8+w0uz5BEEQBN3DyjTzCYIgCLqEMD5BEATBsNOTxkfSaZIekHRrIW0dSZdKuiO/rj2SfaxKg7HsI+k2SX1ZkmhU0GAsX5d0u6SbJZ0rafJI9rEqDcZyVB7HHEmXSHrOSPaxKmVjKTz7hCRLqnqWZURp8L0cIeme/L3MkbT7SPYxSPSk8QF+QP/Ac58Bfmf7+cDv8vvRwA/oP5Zbgb2AK4e9N4PjB/Qfy6XA1ra3Af4GfHa4O9UmP6D/WL5uexvb04ELgP837L1qjx9QEqixWbThLuUHlAedPN729Hw1jAwaDB89aXxsXwnUi5C+BTg9358OvHVYO9UmZWOxPdf2X0eoS23TYCyX2F6W315NCp3R9TQYyyOFtxNoHKW3q2jw/ws0jzbcdQwwlqDL6Enj04ANbN8LkF+njHB/gv68H/jNSHdiMEj6iqR/Ae9i9Mx8+lEx2vBo4pC8JHraaFly73VWJuMTdDGSPg8sA84c6b4MBtuftz2VNI5DRro/7SBpDVK04VFrPOs4EdgCmA7cCxw3st0JYOUyPvdL2gggvz4wwv0JMpIOBN4MvMu9c/DsLGDvke5Em2xBtWjDowLb99tebrsPOBnYYaT7FKxcxud84MB8fyBw3gj2JchI2g34NLCn7cdHuj+DQdLzC2/3BG4fqb4MBtu3VIw2PCqo/dGZeRvJYScYYXpS4UDSj4HXk6TO7we+BPwS+CmwKcl7Z5/REBm1wVgeAr4NrA8sAubYfuNI9bEqDcbyWWA88GDOdrXtD41IB1ugwVh2B14I9JFCe3zI9j0j1ceqlI3F9qmF5/OAGba7PsRCg+/l9aQlNwPzgINr+7/ByNGTxicIgiDoblamZbcgCIKgSwjjEwRBEAw7YXyCIAiCYSeMTxAEQTDshPEJgiAIhp0wPkEQBMGwE8ZnJUPSo0NQ556SPpPv3yrpxW3UcXkr4SFy/r9mDbL6Z9PKwgP0KpLeWwzfIOlMSQ9JevtI9isIBiKMTzBobJ9v+5j89q1Ay8anTd5l+/yhbEDS2KGsv0O8F3jG+Nh+F0nRIwi6ljA+KylKfF3SrZJukfSOnP76PKv4eQ7ydqYk5We757SrJH1L0gU5/b2S/lfSq0iyMl/PQbu2KM5oJK2XT8sjaXVJP8lKw2cDqxf6tqukP0u6QdLPJE2sMJ6XSrpJ0p+BjxbSx+ZxXpfbOjinj5H0nRyU7wJJF9ZmCpLmSfp/kq4C9snjuEjS9ZL+IGnLnG99Sefkuq+TNDOnv64QuOxGSWsO0O9PFvp2ZCH9l7m92yTNKozlB4Xv7PDc5xnAmbm91Ru1FQTdxCoj3YFgxNiLJDmyLUmK5DpJteB02wFbAQuAPwIzJc0Gvge81vZdWcZkBWz/SdL5wAW2fw6Q7VYZHwYet72NpG2AG3L+9YAvADvbfkzSp4GPAV9uMp7vA/9p+wpJXy+kfwBYbPtlksYDf5R0CfBSYBrwElJ4jbnAaYVyT9p+de7T70hSOXdIejnwHWAn4ARSkLKrJG0KXAy8CPgE8FHbf8yG88myDkvaFXg+SehSwPmSXptj0rzf9kPZmFwn6Zzc341tb53LT7a9SNIhwCdsz27yGQVB1xDGZ+Xl1cCPgnF9WgAAA0lJREFUbS8nKX5fAbwMeAS41vZ8AElzSD96jwL/sH1XLv9jYNYg2n8t8C0A2zdLujmnv4K0bPfHbLjGAX8eqCJJk4DJtq/ISWcAb8r3uwLbFPY/JpF+8F8N/CwrHd8n6bK6as/OdU8EXgX8rGBIx+fXnYEXF9LXyrOcPwLflHQm8IvaZ1nCrvm6Mb+fmPt2JfBfkt6W06fm9L8Cz5X0beDXwCUDfS5B0M2E8Vl5aTglAZ4q3C8n/TsZKP9ALOPZ5d3V6p6VCQsKuNT2/i20oQZ11Z79p+2LV0iU/qNJnY/l1zHAohwau54xwCttP1GXfoykX5OERq+WtLPtMoVrAf9t+3t1fXs9ybC90vbjki4HVrP9sKRtgTeSlhb3JQXgC4JRR+z5rLxcCbwj7yOsT5qJXDtA/ttJf3VPy+/f0SDfEqC4xzGPtMQFUPS+upIU7RNJWwPb5PSrSct8z8vP1pD0goEGYnsRsFjSq3PSuwqPLwY+LGnVXN8LJE0ArgL2zns/G5CUj8vqfgS4S9I+ubyyAYA083gmYJyk6fl1ixyW4GvAbGDLBl2/GHh/bU9L0saSppBmZw9nw7MlaTZYW5IcY/sc/n9796/SQBDEcfw7kDKdpaW1jQgWFj6CzyBiYyE2prNKK1gFFSzSiBaKjYVYSEDwTxk1aXwCEQTBxmosdgJHSNQguVz096mOy3HsXnFzuzNhYBOYift0P3ORwlPw+b9OgXugCVwCla/6tcTX/SpwHon4Z+Ctx6VHwEYk2qeALdLL/5qUW+rYAcqx3VYhAp+7v5Cqtw7jt1v6v7yzloBaFBxkVyL7QJvUDO2RlLcqASekPjWdc3d95gMpmC2bWRNoAYtxfg2YjWKBNtBpBbEeRQHNGEvP1uDufkFqOndjZg/AMSmInAOlmH81ngHAJNCIrdA6qR0FcbyrggMZJ2qpID9mZmV3f7eU5KgBT+6+PaKxNPhlkj0znwlS8Jsf14Zp3cysTqbwQ6RotPKRQazEV3eLtDW09831w/QK1K3Hn0wHcBbzuQKqfyjwHAAL9KmyEykCrXxEhszMpkkVeFkf7j43ivGIFIGCj4iI5E7bbiIikjsFHxERyZ2Cj4iI5E7BR0REcvcJgZLjXTzzcEQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ds.tmpprs.isel(time=0,lev=-1).plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's look at a profile"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x20cc84877c8>]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAEWCAYAAABIVsEJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdd3hUddbA8e9JhRBa6JBA6EgvoYi61lVRV+xdwIbu6qq77mtv2Lvr7lpWBbtiQcWCXXRFpYtACL2GmhBKSCD1vH/cGxiGTCZlJjOTnM/zzJOZ2+bcmcmcub8qqooxxhhTVVGhDsAYY0xksgRijDGmWiyBGGOMqRZLIMYYY6rFEogxxphqsQRijDGmWiyBRAAR6Sgie0QkOtSxmMAQkVQRURGJCXUs4UBEfhCRK0Mdh6kaSyBhSETWisgJZY9Vdb2qJqpqSSjj8kVE4kVkooisE5FcEflNREZ5bXO8iCwVkXwRmS4inTzWPSEiK9x9l4rIGK99B4rIPHffeSIy0E88Pp/LXX+CiMwXkTwR2SAi5/k5t0kisltEtojI331sd4+bEE4ob31NiMgxIpIZ6OO6x35VRArdHyhlt2iP9Qki8pyIZIvILhH5XwXHShKRj9zXdZ2IXBSMmGuDiJwnIr+4n6Efyllf6c9kXXpdvFkCMYEQA2wAjgaaAncB74lIKoCItAQ+dJcnAXOBdz32zwP+5O47FnhGREa6+8YBU4E3gebAa8BUd/kh/D2XiPQG3gbucJ9vIDCvgnO7F+gOdAKOBW4WkZO9nrMrcA6wuYLjhLPH3B8oieX8UHkR53U8zP37twqO8yxQCLQBLgaeF5E+wQo6yHKAfwKPeK+o6meSuvW6HExV7RZGN+ANoBTYC+wBbgZSAQVi3G1+AB4AfnG3+RRoAbwF7AbmAKkex+wFfIPzT7EMOK8WzmMhcLZ7fzzwi8e6Ru759fKx7yfATe79E4GNgHisXw+c7GPfCp8LJ3ncX4Xz2Aic6PH4fmCy1zZfAKcAa4ETKnlc7/f0MiADyAVWA1d7xV/qvtd7gPYBfJ9eBR7wsa6n+3lqUonjNML5kuzh9Vl+pJJx/ABc6d6PAu4E1gHbgNeBpl6v21j3c5AN3BHEz/GVwA9eyyr9mazp6xLuN7sCCTOqeinOh/FP6vwafMzHphcAlwIdgK7Ar8ArOL8SM4B7AESkEU7yeBtoDVwIPOfrF5BbXLHTx21hZc5BRNoAPYB0d1Ef4HePc8wDVrnLvfdtCAz12nehuv95roXl7VvJ5xrhPs8iEdksIm+KSJKP82gOtPc8nnu/j8c25wKFqjrNRzyVtQ04DWiCk0yeFpHBbvyjgE164AphUzmx3lrB+7bTz3P/RURy3KKYsz2WD8f5Ep/gFmEt8lrvqQdQoqrLPZYd9FpVwTj3dizQBUgE/uO1zZE4Ce544G4ROay8A9XwdfGlKp/JQL4uYccSSOR6RVVXqeounF/Aq1T1W1UtBt4HBrnbnQasVdVXVLVYVecDU3CKXA6hqn9R1WY+bv39BSUisThXQq+p6lJ3cSKwy2vTXUDjcg7xAs4/2FfV2Lcy2yfjJN6zcYqmGgL/ruBYZfsfciwRSQQeAm70sX+lqern7vupqvoj8DVwVBX2f6SC961ZBbv+C+d1aI1T7PeqiBzhrksG+uKcc3vgOuA1H1/WVX2fKnIx8JSqrlbVPcBtwAVycIODCaq6V1V/x/m8DCjvQDV4XSpSlXMN5OsSdiyBRK6tHvf3lvO47MuvEzDc61fXxUDbQAckIlE4l+eFOF82Zfbg/LL21ASnuMZz/8dxvrDO8/h153NfOdA6bY+I7Knkc+3FSb7L3S+nh3CKnxCRFzyOd7t7LLyO53msCcAbqrqmnJejSkRklIjMdK8Edroxtazpcf1R1fmqut39cTENJ/mf5a7eCxThFHEVuoltOk4RjrdKvceV1B7nyqfMOpx6tjYey7Z43M/nwOe9NlTlXAP5uoQdSyDhKZBDJG8AfvT65ZWoqn8ub2OvL1HvW3p5+7j7CTAR55/8bFUt8lidjscvRLdYrSsHiqkQkQk4RTUnqupur337u8cv0x9I1wOt0xJVNdFj+4qeayE+Xl9VvcbjeA+p6g6cinHPX7cDPI51PHC9OK2ztgApOI0HbvHxMpVLROJxrgqfANq4v4ynAWXn7PfzICK3V/C+7fG3vwf1eN5KFVm6lgMxItLdY5nna1UVm3B++JTpCBRz8I+kSgng6+LJ52eynG0D+bqEn1BXwtjt0BswExjv8TiVQyvRr/RY/wDwqsfjE4CV7v3GOL/gLgVi3dtQ4LAAx/yCG3diOeta4Vy2nw00AB4FZnqsvw1YAbQrZ984N/4bgHicK5t1QJyPOPw91+XAGpyy9QTgPZyrCF/n9QjwI05rm144CeVkd10LnCu5stsG4Nyy1wCnBdcPPo67/z1136MSnFZsgpNI83Ert93n3YtbkRzg9+0cnF/vUThXFrnAMe66WGAlTtFWDHCEu95X44fJwDs4FcdHuO9DH6/zTfWx7/7PNE7F9QqgsxvbB8Cb5f0vlPf/EKDXJdr9/FwD/M+9H1vNz6TP1yXSbyEPwG7lvCkwGqcifSfwD+9/Gu9/GCpIIO7jnsDnQBawHfgeGBjAeDu58e3jQEuhPcDFXjEtdb8If+DgVmIKFHjte7vH+kE4TW33AvOBQX7i8flc7voJ7muRhVPk1ryCY8UDk3BaI20F/l7BtmvxaIWFc0X2oI9tvd/Ta93j73RjmoxH6yg3hu3u+kC2wvrJ/ULbjVOXcIHX+j44DTTygCXAmR7rbge+8HicBHzsbrseuMhj3VHu6xPrI479n2mcZHY3TkLOwm0uW97rVt7/Q4Bel3Hu83jePP/HfH4mq/K6RPpN3BM0xgSYiCwAjlfV7aGOJdRE5E4gS1X/G+pYTOBYAjHGGFMtVolujDGmWiyBGGOMqRZLIMYYY6qlXg0l3bJlS01NTQ11GMYYE1HmzZuXraqtvJfXqwSSmprK3LlzQx2GMcZEFBFZV95yK8IyxhhTLZZAjDHGVIslEGOMMdViCcQYY0y1WAIxxhhTLZZAjDHGVIslEGOMMdVSr/qBGGNMfaCq5OQVsiY7j9VZeazK3sPfTuhBg9jogD6PJRBjjIlwv6zK5rf1O1mdlcfq7D2szspj194Dk4LGRUdxXloKXVsFduZfSyDGGBOhsnILuOeTxUxb5EwR37ZJA7q0asRp/dvRpVUiXVo1omvLRDo0b0h0lPg5WtVZAjHGmAijqkxdsIl7P00nv7CEm0/uyZjDU0mMr92vdEsgxhgTIqpKVm4Bq7Ly2LxrL33aN6VHm0REfF8tbNm1jzs+WsR3S7cxuGMzHjunP91aN67FqA+wBGKMMbVgy659zFu3g9VZe1idnef8zcojt6D4oO1aNY7nyG4tnVv3lrRp0gBwks17czfwwGcZFJWWctdpvRk3MjUoRVOVZQnEGGOCaE9BMc9OX8nEn9ZQWFIKQPumDejaOpGzBnegS6tEurZKpHWTeBZs2MmMFdn8b3kWH/22EYDurRM5sntLVm7bw08rshneOYlHz+5PastGoTwtwBKIMcYERWmpMmV+Jo99tYys3ALOGtSBy47oTNfWjUiIK/+rt0ebxpyXlkJpqZKxZTc/r8zmpxXZvD1rPTFRwv1n9OXiYR2JCuFVhydLIMYYE2Dz1uUw4dMlLMzcxcCUZrx46RAGdWxe6f2jooQ+7ZvSp31Txv+hK/uKSlCFhnGB7cdRU5ZAjDEmQDbt3MsjXyzlk9830aZJPE+fP4DRAzrU+Ioh0B0AA8USiDHG1NDewhL++79VvPDjKlThr8d145qju9KolpvV1ra6fXbGGBNEqsqnCzfzyLQMNu3ax6n92nHrqF6kJCWEOrRaYQnEGGOqYVHmLiZ8ms7cdTvo3a4JT58/kOFdWoQ6rFoVNglERFKA14G2QCnwoqo+IyIDgReABkAx8BdVnS1OT5tngFOAfGCcqs4PTfTGmPpiW+4+Hv9yGR/Mz6RFozgeOasf56alhLQ/RqiETQLBSQ43qep8EWkMzBORb4DHgAmq+oWInOI+PgYYBXR3b8OB592/xhgTcIXFpUz6eQ3//m4FhSWlXHVUF647rhtNGsSGOrSQCZsEoqqbgc3u/VwRyQA6AAo0cTdrCmxy748GXldVBWaKSDMRaecexxhjAmbm6u3c+fFiVm7bwwmHteaOU3vTOQw68oVa2CQQTyKSCgwCZgE3Al+JyBM4E2CNdDfrAGzw2C3TXXZQAhGR8cB4gI4dOwYzbGNMHZO9p4CHpmXw4fyNJDdvyKRxaRzXq02owwobYZdARCQRmALcqKq7ReQB4G+qOkVEzgMmAicA5RU46iELVF8EXgRIS0s7ZL0xxngrLVXembOex75cRn5hMdce25Xrju0edh35Qi2sEoiIxOIkj7dU9UN38VjgBvf++8DL7v1MIMVj92QOFG8ZY0y1pG/axZ0fL+a39TsZ3jmJB8/sG7LRbsNd2CQQt1XVRCBDVZ/yWLUJOBr4ATgOWOEu/wS4TkQm41Se77L6D2NMde0pKOapr5fz6i9raJ4Qx1PnDeDMQR0qHFq9vgubBAIcAVwKLBKRBe6y24GrgGdEJAbYh1ufAUzDacK7EqcZ72W1G64xpi5QVb5YvIUJn6azLbeAC4d15JaTetE0of62rqqssEkgqjqD8us1AIaUs70C1wY1KGNMnbZuex53T03nx+VZ9G7XhBcuqdqgh/Vd2CQQY4ypLQXFJfz3x9U8O30lMVHCXaf1ZuzhnYiJjgp1aBHFEogxpl75ZWU2d368mNXZeZzarx13ndabtk0bhDqsiOQ3gYhIWwBV3SIirYCjgGWqmh7s4IwxJlC27yngwc8z+PC3jXRqkcCrlw3lmJ6tQx1WRKswgYjI1cCtzl15FBgHpAMPi8hjqjox+CEaY0z1qSofzt/IA58vIXdfMdcd243rjusWtnNsRBJ/VyDXAX2AhsA6oJt7JdIcmI7T7NYYY8LSuu153PHRYmaszGZQx2Y8clZ/era1Ph2B4i+BFKlqPpAvIqtUdQuAqu4QEevVbYwJS0Ulpbz80xr++e1yYqOjuH90Hy4e3ils5hKvK/wlkFIRiVXVIuDUsoUi0gBnXCpjjAkrCzbs5NYpC1m6JZeT+rRhwul9rZI8SPwlkLNwx5dS1UyP5S2Am4IVlDHGVNWegmKe/HoZr/6yltaN43nhkiGc3LdtqMOq0ypMIKq63sfyjcDGoERkjDFV9F3GVu76eDGbd+/j0hGd+MdJPev1PB21xfqBGGMi1rbd+5jw6RI+X7SZHm0S+eCikQzpZD3Ja4slEGNMxCktVSbP2cDDX2RQUFzKP07swfg/dCUuxqpma5MlEGNMRFmxNZc7PlrM7LU5jOiSxENn9qNLq8RQh1Uv+etIeLKqfunebwo8BQwFFuNM8rQ1+CEaY4xTSf7Mt8t55ee1NIqP4bFz+nPukGQbbj2E/F2BPAR86d5/Eme62D/htM76L3BG8EIzxhinJ/knv2/iwc8zyNpTwPlpKdx8ci+SGsWFOrR6rypFWGmqOtC9/7SIjA1GQMYYU2bZllzunrqYWWty6J/clBfHpDEwpVmowzIufwmktYj8HWeejiYiIu48HGAdCY0xQbJ7XxH//GYFr/26lsYNYnjozH6cPzSFaOtJHlb8JZCXgLKBY14DWgJZ7gi9C3zuZYwx1aCqfPTbRh6atpTtec7sgP93Yk+aW3FVWPLXkXCCj+VbgDFBicgYUy8t2bSbez5ZzJy1OxiQ0oxJ49Lon2zFVeHMmvEaY0Jq194inv5mOa//upZmCXE8enY/zh2SYgMfRgBLIMaYkCgtVabMz+TRL5eyPa+QS4Z34qYTe9AswYqrIkXYJBARSQFeB9oCpcCLqvqMiLwL9HQ3awbsLGsNJiK3AVcAJcD1qvpV7UdujKmqjM27ueOjRcxfv5PBHZvx6mXD6NuhaajDMlVU7QQiIm3L5gdxHw+uxG5FqrrIx7pi4CZVnS8ijYF5IvKNqp7v8RxPArvc+72BC3AmvGoPfCsiPVS1pJqnZIwJMlXltV/W8tC0pTRuEMPj5/Tn7MHJVlwVoWpyBTIRjzlCgB+BOThNfn3pDKSWt0JVN+N0VERVc0UkA+gALAFnTl3gPOA4d5fRwGRVLQDWiMhKYBjwazXPxxgTRDvzC7n5g4V8vWQrx/VqzRPnDrDOgBGu2glEVU/1WjRHVY8rd2OXiHxfmWOLSCowCJjlsfgoYKuqrnAfdwBmeqzPdJd5H2s8MB6gY8eOlXl6Y0yAzVuXw/XvLGBb7j7uPPUwrjiysw1BUgcErA7EX/Ko7DYikghMAW5U1d0eqy4E3vHctLynKOc5XwReBEhLS7NpeI2pRaWlyvM/ruKpb5bToVlDPrhmJAOsJ3mdUe3e5CLymY/l93k9jhaRtyp5zFic5PGWqn7osTwGZ/ytdz02zwRSPB4nA5sqF70xJtiycgsY+8psHv9qGaP6tuWz64+05FHH1GQ4kqt8LO/oto5CROKBj4AVPrbdz63jmAhkqOpTXqtPAJZ6Tav7CXCBiMSLSGegOzC7iudQKfmFxWzcuTcYhzamTpqxIptRz/zE7DU5PHxWP/594SCbIbAOqnQCEZEkEdk/1Zdb6V2ey4B+bhL5FJiuqvdW4imOAC4FjhORBe7tFHfdBRxcfIWqpgPv4VSyfwlcG6wWWLdMWcQRj3xPfmFxMA5vTJ1RXFLK418t5dJJs2ieEMsn1x3JhcM6Wn1HHeVvPpCOwGPA8cBOZ5E0Ab4HblXVtR7bejbjfQZnuPefgR9FZLCqzq/ouVR1Bj5acKnqOB/LHwQerOi4gZC+aRcACzN3MaJLi2A/nTERaePOvdzwzm/MXbeD89NSuPf0PjSMiw51WCaI/FWivwv8E7i47Ne9iEQD5wKTgREe2z7pte8OoLe7XDnQ/DbiJMY7L1NRSWmIIzEmPH2dvoX/+2AhJaXKMxcMZPTAQxpEmjrIXwJpqaqeFde4iWSyiNzvtfzYQAcXLqLcy++SUmvEZYynguISHp62lFd/WUu/Dk3594WDSG3ZKNRhmVriL4HME5HncIZy3+AuSwHGAr95bigip6lquS2zqrJNOCrrJGsJxJgD1mTncd3b80nftJvLj+jMLaN6Eh9jRVb1ib8EMgZnrKkJOJ30BKf57Cc4LaY8PS4iG6m4J/pDQAQmELsCMcbTjBXZ/PnNeURHCy+NSeOPvduEOiQTAv7mAykEnndv/mwFvJvfevPbnDcclY3TU6qWQIx5b84Gbv9oEd1aJzJx3FA6NGsY6pBMiPhrhRWDcwVyBs4ViOJ01psKTFTVorJtVfWY4IUZWmVFWMV2BWLqMVXlqW+W8+/vV3JU95Y8d/FgGlvfjnrNXxHWGzjNdyfgFF2B0+N7LPAmcL6P/eqUsnmYrQjL1FcFxSXc8sFCPl6wifPTUnjgzL7ERtekH7KpC/wlkMGq2tNrWSYwU0SWBymmsGN1IKY+25VfxPg35jJrTQ7/d1JP/nJMV+sYaAD/CWSHiJwLTFHVUgARicLpB7Ij2MGFC0sgpr7akJPPuFdmsyFnr/XvMIfwdw16AXAOsFVElovICmALzsCGF5S3g4gkiMhdIvKS+7i7iJwWyKBrmzXjNfXRgg07OfO5n8neU8gbVwyz5GEO4a8V1lrceg4RaQGIqmb7OeYrwDzgcPdxJvA+Edh8t8z+KxBrhWXqiS8Xb+HGd3+jVeN43r1sGF1bJYY6JBOG/M4H4s7PcTJOB8Ji9yrk67IirXJ0VdXzReRCAFXdKxFeYBplleimnlBVJv28lgc+X8KA5Ga8PDaNlonxoQ7LhKkKi7BE5DxgOk4CuQ5nythLgQUi0s/HboUi0hB3cicR6QoUBCziELAiLFMflJQqEz5dwv2fLeGk3m1556oRljxMhfxdgdwJjFDVfBFpiTPR00ki0h9ntN2R5exzD87w6inuRFJHAOMCGHOts2a8pq7LLyzm+nd+49uMbVx1VGduG3XY/itvY3zxl0AEKJtJKQ9oDaCqC91h3Q/e2CmqWopTyT7C3f+GStSbhDWxVlimDtuWu48rXp1L+qZd3De6D2MOTw11SCZC+Esg04AvReRHYBROZTgikkQ5Y16pqorIx6o6BPg80MGGilWim7pq+dZcLntlDjl5hbw0Jo3jD7MxrUzl+WuFdYs7K2Bv4D5V/cZdtRMY7GO3mSIyVFXnBDDOkIouqwMpsQRi6o7py7Zx/Tu/0SA2mveuPpx+yU1DHZKJMH5bYanqNJwrEc9lpfiuGD8WuFpE1uEUe4mzi/avYawhs78Vll2BmDpAVXl2+kqe/GY5vdo24aUxQ0hunhDqsEwE8ptAqmFUEI4ZFqwOxES6PQXF/OO93/kyfQujB7bnkbP627SzptqCkUDq3rese0aWQEwkW521h/FvzGNNdh53nnoYVxzZ2ca0MjUSjATyOc5XrgANgM7AMqBPRTuJSArwOtAWKAVeVNVn3HV/xemHUgx8rqo3u8tvwxluvgS4XlW/CsL57M+IlkBMpPouYys3Tl5AbEwUb1w+jJHdWoY6JFMH+JsP5GRV/dK93xRnwqihwGLgb6q61XsfVe3ndYzBwNWViKUYuElV54tIY5zpdL8B2gCjgf6qWiAird3j9sYZj6sP0B74VkR6uHO2B5S6dR+WQEykKS1V/v39Sp7+djl9OzThhUusvsMEjr/BFB/yuP8ksBn4EzAHpyOhX6o6Hyfp+Ntus7stqpoLZOBMYvVn4BFVLXDXbXN3GQ1MVtUCVV0DrMTpKR9wZWnDJpQykSR3XxFXvzmPp79dzlmDOvDBNSMteZiAqkoRVpqqDnTvPy0iY8vbSET+7vEwCqe5b1ZVghKRVGAQMAt4HDhKRB4E9gH/cJsIdwBmeuyW6S7zPtZ4YDxAx44dqxLGfmWNr4pKfA3/ZUx4WbltD+PfmMu67fnc86fejBuZavUdJuD8JZDWbkIQoImIiOr+tqy+rl4ae9wvxqkTmVLZgNzBG6cAN6rqbnda3eY4PduHAu+JSBfK6chIORX4qvoi8CJAWlpatS4h9l+BWD8QEwG+Tt/C39/7nfiYKN66cjgjurQIdUimjvKXQF7iQEJ4DWgJZIlIW2CBj32WqOr7ngvcSane97G953axOMnjLVX90F2cCXzoJq7ZIlLqxpGJM0JwmWSc+doDrixnFpXaFYgJX6Wlyj+/W8G/vlvBgOSmPH/JENo3axjqsEwd5q8n+gQfy7cAY3zsdhuHJovylh3EHUdrIpChqk95rPoYOA74QUR6AHFANvAJ8LaIPIVTid4dmF3Rc1SXXYGYcLdrbxF/f3cB3y3dxrlDkrn/jL40iLX+HSa4KjMfSC+cuoVZqrrHY/n+Flru41HAKUAHEfmXxyGa4BRl+XMEzlDxi0Sk7OrmdmASMElEFgOFwFj3aiRdRN4DlrjHvzYYLbCA/Rmk2K5ATBhasTWX8W/MY0NOPveP7sMlIzpZfYepFf6a8V4PXIvTImqiiNygqlPd1Q/hDNteZhMwFzgdZ0bCMrnA3/wFoqozKL9eA+ASH/s8CDzo79g1pW4GKbIrEBNmpi/bxnVvzadhXAzvjB/B0NSkUIdk6hF/VyBXAUNUdY/bMuoDEUl1O/gd9GWvqr8Dv4vI26paFJRoQ6Ss2UCxtcIyYWTqgo3c9N7v9GrXmJfHDKVt0wahDsnUM/4SSHRZsZWqrhWRY3CSSCd8Xy2kisjDOCP47v9Eq2qXAMQbEgea8doViAkPb8xcx91TFzMsNYmXx6bRuEFsqEMy9ZC/joRbRKSs7wduMjkNpxWUryltXwGex6mXOBZneJI3ah5q6BwowrIrEBNaZSPp3vXxYo7v1ZrXLh9mycOEjL8EMgbY4rlAVYtVdQzwBx/7NFTV7wBR1XWqei9OK6qItb8Iy3qimxBSVR6alsHjXy3jzEEdeP6SIdbSyoSUv2a8mQAi0gqnn0UxsEZV96jqzz522yciUcAKEbkO2Ig7FW6kOtCM165ATGgUl5Ry+0eLeG9uJuNGpnL3ab1tznITcv5aYfUG/gWkAh2B33B6p/+IM9f5rnJ2uxFIAK4H7scpxip32JNIYXUgJpQKiku44Z0FfJm+hRuO786NJ3S3ZromLPgrwpqE07+iG3AksFRVOwM/43T6O4iIRAPnuVcomap6maqeraozvbeNLE7isH4gprblFRRz+atz+DJ9C3ef1pu//bGHJQ8TNvwlkIaqugxAVWfjVpyr6ks4rawO4nbkGyJ17BN+oBmvXYGY2rMjr5CLX57FzNU5PHnuAC4/snOoQzLmIP6a8a4SkbuA74CzcMe/cses8rXvb8BUEXkfZ050ADzGtoo4ZWmj0OpATC3Zunsfl06cxdrt+Tx/8WBO7NM21CEZcwh/CeRynOFEbgd+B25wlyfgeyysJGA7B7e8UiByE4h7CWJXIKY2rNuex8Uvz2JHXiGvXjaUkV1t9kATnvy1wtoJ3FzO8l0cPBeH57rLAhNa+DgwoZRdgZjgyti8mzGTZlNcUso740fQP7lZqEMyxid/dSBVJiI9ROQ7d/BDRKS/iNwZ6OepTdYKy9SGeetyOP+/vxItwvvXHG7Jw4S9gCcQnDlEbgOKAFR1Ic7c5RHL+oGYYPtxeRaXvDybFonxfPDnw+nWurH/nYwJsapMaVtZCao626shVmWGcw9bByaUsisQE3hfLt7CX9+ZT/fWjXnt8mG0ahwf6pCMqRR/HQljgCuAM3EmbVKcYdunAhN9jLqbLSJd3W0RkXOAzYEMOlTsCsQE2uw1OVz/zm/07dCUVy8bRtOGNq6ViRz+rkDeAHYC9+JMIQvOkCZjgTeB88vZ51qcOch7ichGYA1wcSCCDZWyOpBShZJSJdqGkDABsHJbLle9PpfkpIZMGjvUkoeJOP4SyGBV7em1LBOYKSLLy9tBVVcDJ4hIIyBKVXMDEGdIKQeKropKSomOsgHsTM1sy93H2ElziI0WXrtsGM0bxYU6JGOqzF8l+g4ROdcdHBEAEYkSkfOBHeXtICIt3Cltf8KZx/wZEWkRuA89RKgAACAASURBVJBrX5smBybqsRF5TU2VDU+Sk1fIpHFDSUlKCHVIxlSLvwRyAXAOsFVElovICpzh3c/Cd8uqyUAWcLa7bxbwbmDCDY2nzhvI3ac5I7dYPYipieKSUq57ez5LNu3m2YsHWVNdE9H8dSRci1vP4V5FiKpm+zlmkqre7/H4ARE5o0ZRhoHYaKfew/qCmOpSVe6aupjpy7J46Mx+HNerTahDMqZGKt0PRFW3A41F5CwR6VXBptNF5AK3qCtKRM4DPvd3fBFJEZHpIpIhIukicoO7/F4R2SgiC9zbKR773CYiK0VkmYicVNlzqY6YaOelst7oprqenb6Sd2Zv4Npju3LR8I6hDseYGqswgYjIxx73RwPfA38CPhGRcT52uxp4Gyh0b5OBv4tIrojsruDpioGbVPUwYARwrTsfCcDTqjrQvU1z4+mNU4zWBzgZeM4dTj4oYtyWVzYelqmOKfMyeeLr5Zw5qAP/ONG7XYoxkclfK6xOHvdvAY5T1TUi0hJnhN5XvXdQ1Wp1oVXVzbj9RVQ1V0QygA4V7DIamKyqBcAaEVkJDAN+rc7z+xMX4+RaG5HXVNWMFdncMmUhI7u24NGz+9t8HqbO8JdAPH9ux6jqGgBVzRYRn9+kItIfZxbD/cevynDuIpIKDAJmAUcA14nIGGAuzlXKDpzk4jmgYyYVJ5waiYlyEkiRJRBTBRmbd3PNm/Po1jqRFy4dsv+HiDF1gb8EMsAtdhIgXkTaquoWEYkDyi0uEpFJQH8gHSj7tq30cO4ikghMAW5U1d0i8jzO1Ljq/n0SZ5j58n7GHVK+JCLjgfEAHTtWv9y57B+/qNiKsEzlbN61l8temUNifAyvXDaUJg2so6CpW/y1wvJVp5CAU9dRnhGqeshshZXhTlQ1BXir7IpFVbd6rH8J+Mx9mAmkeOyejDPMykFU9UWcnvGkpaVV+9u/rBWWFWGZyti9r4hxk+aQV1DMe9ccTrumDUMdkjEBV+nraRFpLiKNwZknRFV91TX86lH5XWnuNLgTgQxVfcpjeTuPzc4EFrv3PwEuEJF4EekMdAdmV/V5K2t/HUixJRBTscLiUv785jxWZe3hhUuHcFi7JqEOyZig8DeYYnvgEZwK60Rgo1sBOAl40Mdgiq/hJJEtQAFOUZOqan8/sRwBXAosEpEF7rLbgQtFZCBO8dRa3CsfVU0XkfeAJTgtuK5152QPivgYqwMx/qkqt05ZyM8rt/PkuQM4opvNJmjqLn91IG8C96nqGBE5CzgKuBNnvo9ncesWvEzCTQQcqAPxS1VnUH69xrQK9nkQeLCyz1ETsdF2BWL8e/Lr5Xz420Zu+mMPzh6SHOpwjAkqfwmkhar+AE4rKhG5Q1XzgDtFZKmPfdar6ieBDDIcWDNe48/bs9bzn+kruWBoCtcd1y3U4RgTdP4SSJaIXILTgfBsnCKksvoKX/UnS0XkbeBTnCIsoGrNeMNR2RWIFWGZ8ny/dCt3fryIY3q24oEz+lpfD1Mv+EsglwNPALcCC4Dr3OVJOMVY5WmIkzhO9FhW6Wa84SrOTSAFVoRlvCzfmstf3/6N3u2b8OxFg/cPe2NMXeevGe964Lxylm/HaW5b3j6XBSa08BJvrbBMOXblFzH+9bk0jIvh5TFDaRQfjFmijQlPAf+pJCI9ROQ7EVnsPu4vIncG+nlqmxVhGW8lpcr1k39j4869vHDJYNo2beB/J2PqkGBca7+EU7xVBKCqC/E9d0jEsH4gxtsTXy/jx+VZTDi9L2mpSaEOx5haF4wEkqCq3h36ioPwPLUqzvqBGA+f/r6J539YxcXDO9rQ7KbeqlKBrYgciTPi7WJV/drHZtki0hV3XCoROQd3lN1IVjacu12BmCWbdvN/H/xOWqfm3POnPqEOx5iQ8TcfyGyP+1cB/wEaA/eIyK0+drsW+C/QS0Q2AjcC1wQm3NAREeJioiiwK5B6LSevkPFvzKVZwzieu2Swja5r6jV/VyCew4eOB/6oqlki8gTOUOqPlLOPquoJItIIiHLn9ugcoHhDKi46ykbjrcfK5jPfllvA+1cfTuvGVmlu6jd/P5+i3EEUy+ZDzwJwe6P7qteYUraNqua6yz4ISLQhFhcTRWFJ0IbbMmHuoWlL+WXVdh46sx8DUpqFOhxjQs7fFUhTYB7ugIge84Ek4jVulTtPeh+gqTtuVpkmQJ34qRYXHWV1IPXUlHmZTPp5DZcdkco5NsaVMYD/joSpPlaV4gyt7qkncBrQDGfe9DK5wFXVjC+sxMYIRTYner2zMHMnt320iMO7tOD2Uw4LdTjGhA1/w7k3U9Wd3stVNR9Y47VsKjBVRA6vYK6QiGZXIPVPVm4BV78xj1aJ8fznokH7O5QaY/zXgWSLyLcicoWIVKrQt64mD4C4mGgbC6seKSwu5S9vzWNHfiEvjhlCi8T4UIdkTFjxl0AygH8CxwGrRGSqiFwgIvVyfs64aLGOhPXIfZ+lM2ftDh47ZwB92jcNdTjGhB1/CaRIVT9T1Ytx5hx/C2dwxUx3yPZ6JS7GirDqi3dmr+fNmeu5+ugunD6gfajDMSYs+WuFtb+llaruBd4D3hORpsAZ5e4gkoYzc2F7YC/OHObfqmpOQCIOodjoKCvCqgfmrcvh7qmLOap7S24+qVeowzEmbPm7AnmrvIWquktVX/NcJiLjRGQ+zkCKDYFlwDbgSOAbEXlNRCJ60KCGsdHsK7J+IHXZ1t37uObN+bRv1pD/XDiY6CibGMoYX/w1432iCsdqBBzhXqkcQkQGAt2B9VU4ZlhJiI8hv9ASSF21r6iEq9+YR15BMW9dOZymCbH+dzKmHvPXjDcKGIcznW0yTu/zFcALZXOll1HVZys6lqouqEmg4SAhNpr8wogfWNiUQ1W5e+piFmzYyQuXDKZHm8ahDsmYsOevCGsi0BF4GJgOfO4uu1NE/uq5oYjEiMjVIvKliCwUkd9F5AsRuUZE/P6UE5EUEZkuIhkiki4iN3it/4eIqIi0dB+LiPxLRFa6zze4CuddLQnx0eQX2BVIXfTxgo28NzeTvx7XjZP7tgt1OMZEBH+V6EM8pqidISIzVfVuEfkfzhzp//bY9g1gJ3AvkOkuSwbGAm8C5/t5rmLgJlWdLyKNgXki8o2qLhGRFOCPHFz8NQqnSKw7MBx43v0bNI3iYsgvKkFVEbGy8bpiQ04+d32cztDU5tx4Qo9Qh2NMxPCXQIpEpKuqrnJ/4RcCqGqBiHiP6TFYVXt6LcsEZorIcn+BqOpm3HlD3BF8M4AOwBLgaeBmYKrHLqOB11VV3edoJiLt3OMERcO4aEpKlYLiUhrERgfraUwtKi4p5W/vLkCAp88faJXmxlSBvyKs/wOmuwlgivsYEWkFfOa17Q4ROdetN8HdLkpEzgd2VCUoEUkFBgGzROR0YKOq/u61WQdgg8fjTHeZ97HGi8hcEZmblZVVlTAO0SjOSRp7rSK9znjuh1XMXbeDB87sS3LzhFCHY0xE8dcK63sR6QS0UNVsj+VZOFcEni4AHgWeE5GyhNEc+J4qzInujvQ7BWciqmLgDuDE8jYtL+RyzuFF4EWAtLS0Go2EmBDnvFx5hcU0bxRXk0OZMDB//Q6e+W4FZwxsz+iBh/z2MMb44XdKW7eIKLsS263FrefwmD/E736e3Mr2KcBbqvqhiPQDOgO/u3UOycB8ERmGc8WR4rF7MrCpKs9XVQnxzhWINeWNfHsKirlx8gLaNmnAfWf0DXU4xkSkoAwtqqrbVTVbRF6v7D7iZIiJQIaqPuUeZ5GqtlbVVHdo+UycupYtwCfAGLc11ghgVzDrP8CpRAdLIHXBvZ+kk7kjn39eMJAmDay/hzHV4fcKpLJE5BPvRcCxZaP4qurpfg5xBHApsEhEyvqM3K6q03xsPw04BVgJ5AOX+dguYBq6dSD5BdYXJJJ9tnATH8zL5PrjujE0NSnU4RgTsfwmEHfcq5NxKqgVp5joq3LmCUnGaTH1srudAGnAk5UJRFVnUH69huc2qR73Fbi2MscOFLsCiXybdu7l9g8XMSClGX89vnuowzEmolVYhCUiY4D5wDFAAs5wJcfi9NEY47V5Gs70t3fgFCf9AOxV1R9V9ccAxx0SZXUgedYbPSKVlCp/f28BxaXKM+cPtMmhjKkhf1cgd+B0JjzoakNEmgOzgP11HKpaCjwtIu+7f7dW4vgRJcGa8Ua0l35azczVOTx2Tn9SWzYKdTjGRLzKDOdeXtPXUnwUN6lqJnCuiJwK7K5ZeOHlQDNeSyCRZlHmLp78ehmn9GvLuUOSQx2OMXWCvwTyIE6z2a850GmvI86wIvdXtKOqfo4zdladkWCV6BEpv7CYG979jRaN4nnozH42DI0xAVJhIbA750ca8CNQgDOUyQ9Amqq+Guzgwk1sdBRx0VHk25wgEeWBzzNYk53HU+cPoFmCdQA1JlAq05Fwh4hMx6MVlqpWaWiSusQZkdeuQCLF1+lbeHvWeq7+QxdGdm0Z6nCMqVP8zQcyEHgBaIrTiU+AZBHZCfxFVecHP8Tw4swJYlcgkWDb7n3cMmUhfdo34e8n2ii7xgSavyuQV4GrVXWW50K35/crwIDKPImILFLVftWKMMzYrISRobRUuen939lbVMIzFwwkPsZGTzYm0PwlkEbeyQNAVWeKyEHtIEXkLB/HEKBtNeMLO80TYsneUxDqMIwfr/yylp9WZPPAGX3p1tpmFzQmGPwlkC9E5HOc/h5lrbBSgDHAl17bvgu8RfnNfhvUJMhwktw8gdlrckIdhqlAxubdPPrFUk44rDUXD+8Y6nCMqbP8Ded+vYiMwpm8qQPO1UQm8Gw5Y1QtBJ5Q1cXexxGREwIUb8ilNG/I1AV7KSoptZ7MYWhfUQk3Tl5Ak4axPHp2f2uya0wQVaYV1hfAF5U41o347jh4ZlWCCmfJzRMoVdiyax8pSTYBUbh55IulLNuay6uXDaVFYnyowzGmTvM3FlZTEXlERDJEZLt7y3CXNfPcVlV/UtX15R1HVecGMuhQSk5qCDjzaJvwMmNFNq/+spZxI1M5pmfrUIdjTJ3nrwzmPZzpaI9V1Raq2gJnMMWdwPvBDi4cpbjTnmbu2BviSIyn/MJibvtoIV1aNuLWUb1CHY4x9YK/BJKqqo+6EzgBoKpbVPURnCFN6p22TRsQJbBhh12BhJOnv1nOhpy9PHxWPxrEWpNdY2qDvwSyTkRuFpE2ZQtEpI2I3MKBVln1Smx0FO2aNrQrkDDy+4adTJyxhouGd2R4lxahDseYesNfAjkfaAH8KCI5IpKDMxZWEnCe98Yi0ktEjheRRK/lJwco3rCQ3Lyh1YGEiaKSUm6ZspCWifFWdGVMLfM3mOIOVb1FVXupapJ7O8xddlBnCBG5HpgK/BVYLCKjPVY/FPjQQyclKcGuQMLEi/9bzdItudx/Rl+b29yYWlbtjgwi4j0H+VU4k0+dgTOD4V0ickPZ5tV9nnCU3LwhW3P3UVBsQ5qE0uqsPTzz3QpG9W3LSX3qzGAHxkSMmvSEm+D1OFpV9wCo6lqcJDJKRJ6ijiWQlOYJqMJGuwoJmdJS5bYPF9EgJooJo/uEOhxj6iV//UAW+rgtAtp4bb7FHb0XADeZnAa0BPwOpCgiKSIy3e1nkl529SIi97vPuUBEvhaR9u5yEZF/ichKd/3gKp57tXVp5QwDtnzrntp6SuNl8pwNzFqTwx2nHkbrxnVmpBxjIoq/nuhtgJNw+oJ4EuAXr2VjgIMmylDVYmCMiPy3ErEUAzep6nwRaQzME5FvgMdV9S7YX89yN3ANMAro7t6GA8+7f4PusHZNiIkSfs/cycl9reiktm3dvY+Hp2VweJcWnJeWEupwjKm3/CWQz4BEVV3gvUJEfvB87M6FXi5V/dlfIKq6Gdjs3s8VkQygg6ou8disEQcGaxwNvK6qCswUkWYi0s49TlA1iI2mV7vGLMzcGeynMuW4e+piCktKefgsm57WmFDyN5jiFRWsuyjw4ThEJBUYBMxyHz+Ic4WzC6cnPDiDO3r2Rcl0lx2UQERkPDAeoGPHwPV9HJDcjE8WbKK0VImKsi+x2vLl4s18lb6VW0f1IrVlI/87GGOCxl8dSGJF6yu7TVW4x5sC3KiquwFU9Q5VTcEZLv66sk3L2f2QoeRV9UVVTVPVtFatWgUszgEpzcgtKGbN9ryAHdNUbNfeIu6amk6f9k248sjOoQ7HmHrPXyusqSLypIj8wXMCKRHpIiJXiMhXQMA6CYpILE7yeEtVPyxnk7eBs937mThzk5RJBjYFKhZ/BiQ7Y0n+vsGKsWrLI19kkJNXyKNn9yfGhtI3JuT8dSQ8HvgOuBpIF5FdIrIdeBNnlsGxqvpBIAIRpzB7IpChqk95LO/usdnpwFL3/ic4FfTiTrG7qzbqP8p0a51IQlw0CzN31dZT1mu/rtrOO7M3cOWRnenboWmowzHGULn5QKYB3pNHBcMRwKXAIhEpq7S/HbhCRHoCpcA6nBZYuDGdAqwE8gHvjo1BFR0l9E9uypy1NjthsO0rKuG2DxfSqUUCN57QI9ThGGNcfhNIbVHVGZRfr1Fu8nJbX10b1KD8OKp7Kx7/ahnbcvdZX4Qgeua7Fazdns/bVw6nYZyNtGtMuLCC5Bo4uodTKf/T8uwQR1J3pW/axYv/W815acmM7NYy1OEYYzxYAqmB3u2a0DIxnh+XZ4U6lDqpuKSUW6csonlCHLefcliowzHGePGbQETkIvfvBcEPJ7JERQl/6NGSn1ZkUVJ6SAtiU0Ov/LyWRRt3MeH0PjRLiAt1OMYYL5W5AukgIufhNJM1Xo7u0Yod+UUs2mitsQJp3fY8nvxmGX/s3YZT+tlwMcaEI38dCe/BmTzqbSBJRO6ulagiyFHdWyECPy6zYqxAUVVu/2gRsVFR3D+6rw1XYkyY8tcPZAKQA1wC5KjqfbUSVQRJahRH/+Rm/Lh8W6hDqTM+mJfJzyu3c8uoXrRtaq3bjAlXlSnC2qiqk4GNwQ4mUh3doxULNuxkZ35hqEOJeFm5BTzweQbDUpO4aFjgxi4zxgReVVphWS2xD0f3aEWpwk8rrDlvTd37aTp7C0t46Kx+NkilMWHOKtEDYEByU5o2jOUHqwepka/St/D5ws3ccEJ3urUO6BidxpggsEr0AIiJjuKkPm349PdNrMkO/ei8OXmFfJW+hX9+u5xpizaTkxf+RWu78ou48+PF9G7XhPF/6BLqcIwxleBvPpAJIvJ/OJXoyar6RO2EFXn+cWJPvli0hbunLub1y4fVasuhjTv3MmdNDrPX5jB7TQ4rtx061e5h7ZowsmsLRnZtwbDOSTRuEFtr8VXGQ9OckXZfGTeUWBtp15iIUJmxsDap6mQRuTDo0USw1k0acNOJPbj30yV8vmgzp/VvH5TnUVVWZe1h9podzHETxsadewFoHB9DWmpzzhrcgWGpSfRu34SlW3L5ddV2flmVzZsz1zFxxpr9A0E6CaUlQzo1p0Fs6MaYmrEim3fnbuDPx3S1kXaNiSDijElYP6SlpencuXODdvySUmX0szPYtruA7246OiC/8otLSlmyeTez1+QwZ20Oc9fuYLtbJNUyMZ5hnZszLDWJoZ2T6NW2CdEVVDzvKyph/vodbkLZzoINOykpVeKioxjcqRkju7ZkZNcWDEhpVmtXAfmFxZz49P+Ii45i2g1HhTSRGWPKJyLzVDXtkOWWQAJrwYadnPncz4we0J4Jp/elaULVk8jGnXv5Ydk2pi/NYubq7ewpKAagY1ICQ1OTnKTRuQWpLRJqVFS2p6CYOWty+GVVNr+s2s6SzbtRhYS4aIZ1TuKvx3VnSKfm1T5+ZUz4NJ1Xfl7L+9ccztDUpKA+lzGmeiyBUDsJBODxr5by7PRVJMbHcMmITlw8vCPJzRv6/LIvKill7todTtJYto3lW506jOTmDTm6RyuGd2nBsNSkoHeq25FXyMzVztXJ10u2kFdQwuTxI4JWrDRv3Q7OeeEXLh3RiftG9w3Kcxhjas4SCLWXQAAyNu/m2ekr+XzRZlShYWw0nVok0LllI1JbNqJzi0aUqPLjsixmrMxmT0ExsdHCsM5JHNuzNcf0bEXXVokhG8Zj8669nP3cLxSWKFP+fDidWjTyv1MVFBSXcOq/ZrC3sISv/vYHEuPDZmoaY4wXSyDUbgIpsyY7jxkrsliTnc/a7Xmszc5jfU4+xe7ove2aNuCYnq05tmcrRnZrGVZfpCu35XLOC7/StGEsH1wzklaN4wN27Ce/Xsa/v1/Jq5cN5ZierQN2XGNM4PlKIOHzbVVHdW7ZiM4tD/71XlxSysadeykqUbq2ahS2gwV2a92YV8YN5aKXZjHuldlMHj8iIA0DlmzazfM/rOKswR0seRgTwazBfQjEREfRqUUjurUOXRFVZQ3q2JznLhnM0i25XP3GPAqKS2p0vOKSUm6ZspBmCbHcdWrvAEVpjAkFSyDGr2N7tuaxs/vzy6rt/P3d32s0edbLM9awaOMu7hvdl+aNbJIoYyJZ2CQQEUkRkekikiEi6SJyg7v8cRFZKiILReQjEWnmsc9tIrJSRJaJyEmhi77uO3tIMref0ovPF21mwqfpVKfubHXWHp7+Zjkn9WnDqL42SZQxkS5sEghQDNykqocBI4BrRaQ38A3QV1X7A8uB2wDcdRcAfYCTgedExHqhBdH4P3Rl/B+68Pqv6/jP9yurtG9pqXLrlEXEx9gkUcbUFWGTQFR1s6rOd+/nAhlAB1X9WlWL3c1mcmBU4NHAZFUtUNU1wEpgWG3HXd/cenIvzhrUgSe/Wc47s9dXer+3Zq9n9toc7jytN62b2CRRxtQFYZNAPIlIKjAImOW16nLgC/d+B2CDx7pMd5n3scaLyFwRmZuVZcOt11RUlPDoOf05pmcr7vhoEV+lb/G7z8ade3lkWgZHdmvJuUNsVgBj6oqwSyAikghMAW5U1d0ey+/AKeZ6q2xRObsfUjCvqi+qapqqprVq1SoYIdc7sdFRPHfxYPonN+Ov7/zGrNXbfW6rqtzx0SJKFR4+q58VXRlTh4RVAhGRWJzk8ZaqfuixfCxwGnCxHqi9zQRSPHZPBjbVVqz1XUJcDK+MG0pK84Zc+fpcMjbvLne7jxds5IdlWdx8ck9SkhJqOUpjTDCFTQIR56fpRCBDVZ/yWH4ycAtwuqrme+zyCXCBiMSLSGegOzC7NmOu75o3iuP1K4bTKC6GsZNmsyEn/6D1WbkFTPh0CYM7NmPM4amhCdIYEzRhk0CAI4BLgeNEZIF7OwX4D9AY+MZd9gKAqqYD7wFLgC+Ba1W1Zr3cTJV1aNaQ1y4fxr6iEsZOms32PQX71937aTr5BSU8dk7/CoeZN8ZEJhsLywTE3LU5XPzyLHq2bcw7V41gxspsrn5jHv84sQfXHdc91OEZY2rAxsIyQZWWmsSzFw3m6jfnMf6NuazYuodebRtz9dFdQx2aMSZIwqkIy0S4E3q34eEz+/Hzyu1szyvk8XMG2PzmxtRhdgViAuq8oSn76zv6Jdv85sbUZZZATMCdbZ0FjakXrHzBGGNMtVgCMcYYUy2WQIwxxlSLJRBjjDHVYgnEGGNMtVgCMcYYUy2WQIwxxlSLJRBjjDHVUq8GUxSRLGBdNXdvCWQHMJxIYOdcP9S3c65v5ws1P+dOqnrIjHz1KoHUhIjMLW80yrrMzrl+qG/nXN/OF4J3zlaEZYwxplosgRhjjKkWSyCV92KoAwgBO+f6ob6dc307XwjSOVsdiDHGmGqxKxBjjDHVYgnEGGNMtVgC8UNEThaRZSKyUkRuDXU8gSIiKSIyXUQyRCRdRG5wlyeJyDcissL929xdLiLyL/d1WCgig0N7BtUnItEi8puIfOY+7iwis9xzfldE4tzl8e7jle761FDGXV0i0kxEPhCRpe77fXhdf59F5G/u53qxiLwjIg3q2vssIpNEZJuILPZYVuX3VUTGutuvEJGxVYnBEkgFRCQaeBYYBfQGLhSR3qGNKmCKgZtU9TBgBHCte263At+panfgO/cxOK9Bd/c2Hni+9kMOmBuADI/HjwJPu+e8A7jCXX4FsENVuwFPu9tFomeAL1W1FzAA59zr7PssIh2A64E0Ve0LRAMXUPfe51eBk72WVel9FZEk4B5gODAMuKcs6VSKqtrNxw04HPjK4/FtwG2hjitI5zoV+COwDGjnLmsHLHPv/xe40GP7/dtF0g1Idv+xjgM+AwSnh26M93sOfAUc7t6PcbeTUJ9DFc+3CbDGO+66/D4DHYANQJL7vn0GnFQX32cgFVhc3fcVuBD4r8fyg7bzd7MrkIqVfRDLZLrL6hT3kn0QMAtoo6qbAdy/rd3N6spr8U/gZqDUfdwC2Kmqxe5jz/Paf87u+l3u9pGkC5AFvOIW270sIo2ow++zqm4EngDWA5tx3rd51O33uUxV39cavd+WQCom5SyrU+2eRSQRmALcqKq7K9q0nGUR9VqIyGnANlWd57m4nE21EusiRQwwGHheVQcBeRwo1ihPxJ+zWwQzGugMtAca4RTheKtL77M/vs6xRuduCaRimUCKx+NkYFOIYgk4EYnFSR5vqeqH7uKtItLOXd8O2OYurwuvxRHA6SKyFpiMU4z1T6CZiMS423ie1/5zdtc3BXJqM+AAyAQyVXWW+/gDnIRSl9/nE4A1qpqlqkXAh8BI6vb7XKaq72uN3m9LIBWbA3R3W2/E4VTEfRLimAJCRASYCGSo6lMeqz4BylpijMWpGylbPsZtzTEC2FV2qRwpVPU2VU1W1VSc9/J7Vb0YmA6c427mfc5lr8U57vYR9ctUVbcAG0Skp7voeGAJdfh9xim6GiEiCe7nvOyc6+z77KGq7+tXwIki0ty9cjvRXVY5oa4ECvcbcAqwHFgF3BHqeAJ4XkfiXKouBBa4t1Nwyn6/A1a4f5Pc7QWnhbgVNwAAA89JREFURdoqYBFOC5eQn0cNzv8Y4DP3fhdgNrASeB+Id5c3cB+vdNd3CXXc1TzXgcBc973+GGhe199nYAKwFFgMvAHE17X3GXgHp46nCOdK4orqvK/A5e65rwQuq0oMNpSJMcaYarEiLGOMMdViCcQYY0y1WAIxxhhTLZZAjDHGVIslEGOMMdViCcQYY0y1WAIx9ZKIrBWRlu5Q53/xWN5eRD5w7w8UkVOqcex7ReQfVdx+o4jcV8XneblsdOiy83Hv73H/ep7LOBH5T1WOX87zvSUiOSJyjv+tTX1gCcTUd82A/QlEVTepatkX5ECczpW14WlVvbsqO6jqlaq6pIL1nudSIyISrf/f3v28RBWFYRz/PmGBlATuBcFFEUL2g34solq0iiCLcNGuQAqiZQRBiW0iW7doVYsWUWYEkVSQZFKRiSUtokX9A7ZJBsv0bfEeaxxmHL1Gt8b3s5lh7rn3njkwHM49zPP6v/ZrIokh/BkxgYSaJumepDepuFBnmSaXgBZJo5J6JDWnIkSrgG6gIx3rKF1ZpHbN6f05eeGxJ8C6ojYtkvpTHwYlrV9An7sk3ZD0KK0sDkm6LGksXWtlajcgaes812lWUbEhoCmd/0HShWpjJGlCUrekV3j8eQhz1FVvEsJ/7ZiZfZFUD7yW1Gtm40XHzwKtZtYGv6LtMbPvks7jkQ+n0rGucjeQtAXP1tqE/6ZG8PhwgGvACTP7KGk7cBUPcaymBdiLFzJ7ARw2szOS+oD9eCTJYm0DWoECPhYPzGyYymO0Gq81saiVUVg+YgIJte60pPb0vgmvyDY+T/ssdgF9ZlYAkHQ/va7BU2Bve6Yf4JlMC/HQzKYkjeEV9frT52N4EaEsHs9OnpLu4nlow1Qeo2k8rTmEsmICCTVL0h482nunmRUkDeDBeVn9YO5j3+JrlQuVW4EXMWrLcK9vAGY2I2nKfofWzZD9d1vaR6syRpNmNp3xXmEZiD2QUMvW4rWuC2nvYUeZNl+Bhgrnlx77jNfSQNJmvGARwDOgXVK9pAbgAIB5ga5Pko6kcyRp49K+0pLsk9SYHlUdBIZY2BiFUFZMIKGW9QN1kt4BF4GXpQ3SI52htCHeU3L4KbBhdhMdf5zTKGkUOInH/GNmI8AtPBK/FxgsusZR4Likt8B7vFJeXp7j0eajQG/a/6g6RiFUEnHuIeQsbc5PmNmVvPtSjaTreB2VO3n3JeQvViAh5G8C6FzsHwn/Nkk3gd3AZN59Cf+GWIGEEELIJFYgIYQQMokJJIQQQiYxgYQQQsgkJpAQQgiZ/ARUHNDejSk/wgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ds.tmpprs.sel(lon=10,lat=56).isel(time=0).plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Convert to dfs3\n",
    "\n",
    "Dfs3 does not support irregularly spaced spatial axis as is used by the vertical coordinate axis in this case (pressure levels).\n",
    "\n",
    "**Thus, please note that the vertical coordinates are not correct in this example.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(10.0, 54.0, 21, 17, 0.24, 0.24)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lat = ds.lat.values\n",
    "lon = ds.lon.values\n",
    "lev = ds.lev.values\n",
    "\n",
    "nx = len(lon)\n",
    "ny = len(lat)\n",
    "nz = len(lev)\n",
    "\n",
    "x0 = lon[0]\n",
    "y0 = lat[0]\n",
    "\n",
    "dx = np.round((lon[-1] - lon[0]) / nx,2)\n",
    "dy = np.round((lat[-1] - lat[0]) / ny,2)\n",
    "\n",
    "x0, y0, nx, ny, dx, dy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-06-04T00:00:00.000000000\n"
     ]
    }
   ],
   "source": [
    "t = ds.time.values\n",
    "print(t[0])\n",
    "start_time = pd.to_datetime(t).to_pydatetime()[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Variable types"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Temperature"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from mikeio.eum import EUMType\n",
    "EUMType.Temperature"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[degree Celsius, degree Fahrenheit, degree Kelvin]"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EUMType.Temperature.units"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Data manipulation\n",
    "Flip upside / down"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "temperature = ds.tmpprs.values\n",
    "\n",
    "temperature = temperature[:,::-1,:,:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "from mikeio.eum import ItemInfo, EUMUnit\n",
    "\n",
    "d = []\n",
    "d.append(temperature)\n",
    "      \n",
    "dfs = Dfs3()\n",
    "dfsfilename = f\"gfs_{dtstr}_{hour}_temperature_profile.dfs3\"\n",
    "dfs.create(filename=dfsfilename,\n",
    "           data=d,\n",
    "           start_time = start_time,\n",
    "           dt=3600,\n",
    "           items=[ItemInfo(EUMType.Temperature, EUMUnit.degree_Kelvin)],\n",
    "           coordinate=['LONG/LAT', x0, y0, 0],\n",
    "           length_x=dx, length_y=dy\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clean up (don't run this)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "os.remove(dfsfilename)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
